Hello, this is to address a note on recurrence relations and why there have been two cases now where I computed an “i” term rather than an “i+1” term.

Here I will repeat the Projectile with Gravity and Drag and the 2-Body problem recurrence relationship

Projectile with Gravity and Drag Forces:

2D-2 Body Problem:

For the projectile, p_i and q_i are the terms computed with just “i” and that goes for s2.s3, and a in the in the 2D-2 Body problem.

When setting these equations up, you cannot multiply more than 2 variables at a time. The cauchy product is designed to multiply 2 variables at a time. If there are more(in both cases above there are), another variable is there as a placeholder. In the 2D 2-Body problem, “s” is computed by multiplying “a” and what would be “s^3.” To do this, an “s2” variable is made that is a cauchy product of s and s. Then an “s3” variable is made that is a cauchy product of “s2” and “s” to update “s” completing the recurrence setup. Also, the variable”a” is multiplying a difference, which can be done. The velocity variables are doing that too.

I’ll demonstrate an example that shows this idea.

y’ = y^3

This has a known analytical solution of 1/sqrt((1-2t)). There is a specified domain so that it doesn’t blow up.

A variable needs to be created that is a cauchy product of two y terms so that the final way to update y is multiplying two variables.

let a = y^2

so

y’ = a*y

Now here is the recurrence relationship:

We start with a y(0) that is the given initial condition. Then that is used to compute a(0). This cycle repeats and then you use that to get the Maclaurin Series for only y(t). This is a placeholder and so the initial condition for the next frame or calculation is not the Maclaurin Series.

Below is Matlab Code:

————————————————————-

% Adam Wermus

% November 16, 2015

% This program solves the ODE y’ = y^3 with the analytical solution

% y = 1/sqrt((1-2t)) ; y(0) = 1

% and the Parker Sochacki version

% Clear everything

clc;

close all;

clear all;

% Analytical solution

t = -10:.01:-5; % must be same length as frames

z = 1./sqrt(1-(2.*t));

% Parker-Sochacki Setup

tau = t(2)-t(1); % time step

order = 100; % order of Maclaurin Series

frames = length(t); % total number of frames

y_sum = zeros(1,frames);

a_sum = zeros(1,frames);

a = zeros(1,order);

y = zeros(1,order);

% Parker-Sochacki Implementation

y_sum(1) = 0.2182;

y(1) = y_sum(1);

a_sum(1) = y_sum(1).*y_sum(1);

a(1) = a_sum(1);

% Power Array for time

power_array = tau.^(0:order-1);

for j = 2:frames

% Cauchy product setup

cauchy_product_y = 0;

cauchy_product_a = 0;

% Compute Maclaurin Coefficients

for i = 1:order

cauchy_product_y = 0;

cauchy_product_a = 0;

for k = 1:i

cauchy_product_y = cauchy_product_y + (a(k)*y(i+1-k))/i;

cauchy_product_a = cauchy_product_a + (y(k)*y(i+1-k));

end

y(i+1) = cauchy_product_y;

a(i) = cauchy_product_a;

end

% stores the values of the maclaurin series

y_temp_sum = 0;

a_temp_sum = 0;

% Compute Power series

for i = 1:order

y_temp_sum = y_temp_sum + (y(i)*tau^(i-1));

a_temp_sum = a_temp_sum + (a(i)*tau^(i-1));

end

y_sum(j) = y_temp_sum;

a_sum(j) = a_temp_sum;

% reset initial conditions

y(1) = y_temp_sum;

a(1) = a_temp_sum;

end

% Plot Exact Solution

subplot(2,1,1);

plot(t,z,’r’);

title(‘dy/dt = y^3’)

xlabel(‘t’)

ylabel(‘y’)

legend(‘Analytical Solution’);

% Plot Parker-Sochacki Solution

subplot(2,1,2);

plot(t,y_sum,’r’);

title(‘dy/dt = y^3’)

xlabel(‘t’)

ylabel(‘y’)

legend(‘Approximate Solution’);

% plot(1:frames,zsum); PS method

————————————————————-

Here is a plot with the analytical solution on top and the Parker-Sochacki solution on the bottom.

I hope that clears up why some recurrence relations have only an i term in them.