Subscript indices must be real positive or logicals

      %Astro Final: Problem 2b Competition of Proton Capture and Beta
  %Decay
  %Carrie Elliott - 2015

T = linspace(0.1, 1, 100);
                        %Establishes range(T) in units of 10^9 Kelvin

p1B = -4.538;           %only parameter for Beta decay

Rbeta = exp(p1B);        %The rate of Beta decay is constant, and has only one 
                        %parameter independent of temperature.

p1_p1 = 17.29;          %These are the seven parameters
p2_p1 = 4.341*10^-4;    %for the first component of proton capture;
p3_p1 = -18.10;         % _p1 denotes this.
p4_p1 = 0.1471;         %
p5_p1 = -0.1550;        %
p6_p1 = 0.0167;         %
p7_p1 = -0.7411;        %

Rp1 = exp(p1_p1 + (p2_p1\T) + (p3_p1\((T).^(1\3))) + p4_p1(T.^1\3) + p5_p1*T + p6_p1*(T.^(5\3)) + p7_p1(ln(T)));

p1_p2 = 3.428;          %These are the seven parameters
p2_p2 = -4.654;         % for the second component of proton capture;
p3_p2 = -3.999;         % _p2 denotes this.    
p4_p2 = 8.947;          %       
p5_p2 = -0.6332;        %
p6_p2 = 0.0402;         %
p7_p2 = -5.113;         %

Rp2 = exp(p1_p2 + (p2_p2/T) + (p3_p2/((T)^(1/3))) + p4_p2(T^1/3) + p5_p2*T + p6_p2*(T^(5/3)) + p7_p2(ln(T)));

Rproton = Rp1 + Rp2;

hold on
plot(Rbeta)
plot(Rproton)
title('Competition between Beta+ Decay and Proton Capture as Function of Temperature');
xlabel('Temp, 10^9 Kelvin');
ylabel('Rate');
hold off %

What's going on here? I'm simply trying to plot two rates of reactions over a range of temperatures (0.1*10^9 - 1.0*10^9) on the same plot, and it says the error is in line 21 but I have no idea where or why.

Answer 1

[Star Strider]

This is probably the first instance that throws the error:

p4_p1(T.^1\3)

but p4_p1 is a scalar constant, not an array:

p4_p1 = 0.1471;         %

What do you want to do in that line? My guess is that you want to multiply where you specified no operator, so see if substituting this line solves the problem:

Rp1 = exp(p1_p1 + (p2_p1\T) + (p3_p1\((T).^(1\3))) + p4_p1*(T.^1\3) + p5_p1*T + p6_p1*(T.^(5\3)) + p7_p1*(log(T)));

and:

Rp2 = exp(p1_p2 + (p2_p2./T) + (p3_p2/((T).^(1/3))) + p4_p2(T.^1/3) + p5_p2*T + p6_p2*(T.^(5/3)) + p7_p2*(log(T)));

Also, the natural logarithm function is log not ln, and the two division operators ‘/’ and ‘\’ are not the same, but you seem to have used one in Rp1 in the same place you used the other in Rp2. I did not change the division operators because I don’t know which one you want.

[Carrie Elliott]

@Star Strider Thank you. I fixed all those things and am now getting this:

Error using  + 
Matrix dimensions must agree.

Error in BetaProton (line 21)
Rp1 = exp(p1_p1 + (p2_p1\T) + (p3_p1\(T.^(1\3))) + p4_p1*(T.^1\3) + p5_p1*T + p6_p1*(T.^(5\3)) + p7_p1*(log(T)));

[Star Strider]

My pleasure.

When I took that expression apart, it turns out that this term is causing the problem:

p4_p1*(T.^1\3)

It produces a (100x1) (column) vector, however the other terms are (1x100) (row) vectors.

Since I don’t know what you want to do, I’ll let you sort that. You will probably have to do the same change to ‘Rp2’.

I’ll be here to help with any other problems.

[Carrie Elliott]

Actually, the divides for the exponents (1/3 and 5/3) should've had forward slashes. So thank you, that ultimately fixed it

[Star Strider]

My pleasure.

I thought that might be the problem. I’m glad you got it sorted.