I wrote a code in MATLAB to solve a system of differential equations, but my solution doesn't seem to take into consideration the initial conditions I specified. I am not sure how to interpret this problem. In the code below, r_0 is equal to 32612336, but the output rsol starts at 0: output I will be extremely grateful for any ideas/help!

% %function definitions
syms r(T) p(T) T Y t

G= 4.30091252525*(10^(-3)); %grav constant, in (parsec*km^2)/(Ms*sec^2) 
M = 170000; %mass of the Black Hole, in solar masses
m = 30; %mass of the object, in solar masses
reduced = m*M / (m+M);
c = 0.0020053761; % AU/sec, speed of 
rs = 2*G*M / c^2; %Schwarzchild radius
PI = 3.14;
E = 89656170000.0000;
L = 2.592191000000000e-13;

ode1 = diff(r,T) == sqrt((E^2/c^2 - c^2) + 2*parsec_to_AU(G)*M/r*(1+L^2/((c^2)*(r^2)))-(L^2)/(r^2));
ode2 = diff(p,T) == (1/r^2) * L;

% %function setup
[f1, Subs]= odeToVectorField([ode1 ode2]);
F1= matlabFunction(f1, 'Vars',{T,Y})

r_0= 32612336; %In KM
p_0 = PI;
span = [0 10^(-11)];
conditions = [r_0;p_0]; 
options = odeset('MaxStep', 10^(-14),'AbsTol',10^(-7));
[rsol,psol] = ode45(@(T,Y)F1(T,Y),span,conditions,options);

You're printing out the time steps.

If you print out the solution, you should find the IC is initialized properly.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.