# Index exceeds the number of array elements (5).

4 views (last 30 days)
Hi, I need help how to solve this issue at line67 if(Af(i)<0). When I change for i =[1,2,3,4,5,6] to for i=1:m-1, the problems still persist. Please enlighten me.
clc;
close all;
clear all;
% Constants
rGEO = 42164e3; % [m] Radius of the GEO orbit
rE = 6378e3; % [m] Radius of Eearth
mu = 3.986e14; % {m^3/kg] Gravitational parameter of Earth
sigma = 5.66961e-8; % [W/m^2/K^4] Stefan-Boltzmann constant
S = 1372; % [W/m^2] Solar constant
% Material parameters
c_AL = 896; % [J/kg/K] Heat capacity of Aluminium
r_AL = 2.7e3; % [kg/m^3] Density of Aluminium
a_AL = 0.15; % [-] Absorption of Aluminium
e_AL = 0.05; % [-] Emissivity of Aluminium
k_AL = 229; % [W/m/K] Heat transfer of Aluminium
% Satellite configuration
e_ik = 0.2; % [-] View factor
l = 2; % [m] Length of plate
d = 0.02; % [m] Thickness of plate
mk = l^2*d*r_AL; % [kg] Mass of plate
% GEO Orbit
omega_s = sqrt(mu/rGEO)/rGEO; % [rad/s] Orbital angular rotation speed of the satellite
% Simulation parameters
dt = 60; % [s] Step size
n = 4; % [-] Number of orbits
tE = n*2*pi/omega_s; % [s] Duration of simulation for 90days
% -------------------------------------------------------------------------
%% Initialization
phi = 0.1; % [rad] Anomaly of the satellite
t = 0; % [s] Start time
m=1; % Counter [-]
Tsp(m,:)=273*ones(1,6); % [K] Initial temperature for simulation in spring
Tsu(m,:)=273*ones(1,6); % [K] Initial temperature for simulation in summer
% -------------------------------------------------------------------------
%%Eclipse period
alpha=2*asin(rE/rGEO);
% -------------------------------------------------------------------------
%% Main loop
while t<tE
%% Calculation of temperature with the Euler met
m=m+1;
%%The projected areas for the nodes with respect to the incident Sunlight are
Af(1)=l^2*cos(phi); %Spring node 1
Af(2)=0; %Spring node 2
Af(3)=l^2*sin(phi); %Spring node 3
Af(4)=-l^2*sin(phi); %Spring node 4
Af(5)=0; %Spring node 5
As(6)=-l^2*cos(phi); %Spring node 6
As(1)=l^2*cos(phi)*cos(theta); %Summer node 1
As(2)=l^2*sin(theta); %Summer node 2
As(3)=l^2*sin(phi)*cos(theta); %Summer node 3
As(4)=-l^2*sin(phi)*cos(theta); %Summer node 4
As(5)=0; %Summer node 5
As(6)=-l^2*cos(phi)*cos(theta); %Summer node 6
%Applying Loop over all nodes
for i =[1,2,3,4,5,6]
%For negative areas both spring and summer, we need to make them 0 instead
%of absolute values
if(Af(i)<0)
Af(i)=0;
end
if(As(i)<0)
As(i)=0;
end
%Emitted energy from Sun in summer
Qs_sun = S * a_AL * As(i);
%Emitted Energy from sun in spring
if (phi>(pi-alpha/2) && phi <(pi+alpha/2))
Qf_sun=0;
else
Qf_sun = S * a_AL*Af(i);
end
%Heat emission is conducted on both sides
%Heat emission from all plates are proportional to T to the power of 4
Qf_em = -e_AL * sigma * 2 * l^2 * Tsp(m-1,i)^4; %[W] Heat emission during spring
Qs_em = -e_AL * sigma * 2 * l^2 * Tsu(m-1,i)^4; %[W] Heat emission during summer
Qf_ab =0.0;
Qs_ab =0.0;
Qf_HT =0.0;
Qs_HT =0.0;
for k=[1,2,3,4,5,6]
%Summation of the absorbtion from all other plates
if (i ~= k)
Qf_ab = + Qf_ab + a_AL * e_ik * (e_AL * sigma * Tsp(m-1,k)^4) * l^2;
Qs_ab = + Qs_ab + a_AL * e_ik * (e_AL * sigma * Tsu(m-1,k)^4) * l^2;
end
%Summation of the heat conduction
if (k + i ~= 7)
Qf_HT = Qf_HT - k_AL * d * (Tsp(m-1,i) - Tsp(m-1,k));
Qs_HT = Qs_HT - k_AL * d * (Tsu(m-1,i) - Tsu(m-1,k));
end
end
%Total temperature of spring
Tsp(m,i) = Tsp(m-1,i) + (Qf_em+Qf_ab+Qf_HT+Qf_sun)*dt/(mk*c_AL);
%Total temperature of summer
Tsu(m,i) = Tsu(m-1,i) + (Qs_em+Qs_ab+Qs_HT+Qs_sun)*dt/(mk*c_AL);
end
%There are two loops of variables i and t. It is needed to end them
%respectively. First loop is i loop
phi = phi + omega_s * dt;
if phi > 2*pi
phi = phi - 2*pi;
end
t = t + dt;
end %while
% -------------------------------------------------------------------------
%% Plot of temperature curves
figure('Position',[100 100 1000 800],'Name','Node-model');
subplot(2,1,1)
hold on;
plot(Tsp(:,1),'k');
plot(Tsp(:,2),'b');
plot(Tsp(:,3),'y');
grey=plot(Tsp(:,4)); set(grey,'Color',[0.5 0.5 0.5]);
plot(Tsp(:,5),'g');
plot(Tsp(:,6),'r');
title('Temperature Profile of 6 nodes on Spring');xlabel('Time [min]');ylabel('Temperature [K]')
grid on; legend('Node1','Node2','node3','Node4','Node5','Node6');
subplot(2,1,2)
hold on;
plot(Tsu(:,1),'k');
plot(Tsu(:,2),'b');
plot(Tsu(:,3),'y');
grey=plot(Tsu(:,4)); set(grey,'Color',[0.5 0.5 0.5]);
plot(Tsu(:,5),'g');
plot(Tsu(:,6),'r');
title('Temperature Profile of 6 nodes on Summer');xlabel('Time [min]');ylabel('Temperature [K]')
grid on;legend('Node1','Node2','Node3','Node4','Node5','Node6');

Tommy on 17 Jun 2020
Likely a typo in line 56:
As(6)=-l^2*cos(phi); %Spring node 6
% ^ should be Af?
As it is, Af only has 5 values.
Yes thank you! now solved!