Division by zero error

I have the following functions in which I am doing a Ritz Method analysis for a tapered cantilever beam, under clamped free conditions.

I am repeatedly getting a zero division error even though I have checked that there are no zeros in the denominators.

Can I know how do I fix this?

I have attached the file code and pasted the lines for reference.

Thanks in advance

>> %Constants
E = 2*10^12;                                                        %Young's Modulus
rho = 300;                                                          %Density
    
L = 29;                                                             %Length

syms i j y                                                          %Declare i, j, y variables

%Variable Equations
o = y/L;                                                            %Position Ratios

A = 4.1*(1.005-o);                                                   %Area
I = 0.03*(1.005-0.5*o-0.5*o^2);                                      %Moment of Inertia

%Miscellenous
n = 7;                                                               %Iterations

%Basis Functions Cosine
phis_i = cos((i-3)*pi*o);
phis_j = cos((j-3)*pi*o);

%For Mc Matrix
sig_fun = rho*A*phis_i*phis_j;                                              %Function for integral a
sigma = int (sig_fun, y, 0, L);                                             %Integral a for y from 0 to L

%Matrices Cosine Series
Mc = zeros (n);

for si = 1:n
    for sj = 1:n
        Mc(si,sj) = subs (sigma,{i,j},{si,sj});
    end
end

[KALYAN ACHARJYA]

MATLAB code?

Answer 1

[Walter Roberson]

        %Mc(si,sj) = subs (sigma,{i,j},{si,sj});
        t1 = limit(sigma, i, si);
        t2 = limit(t1, j, sj);
        Mc(si, sj) = t2

[Elton Sim]

Hi Walter, thanks for this. It worked!

Answer 2

[Ameer Hamza]

Division by zero occurs in your equation. Following image show expression for sigma from the live editor

As you can see, there are several conditions when 0 can come in the denominator. One such condition is i==j. Therefore, when the loop starts, you have i=1, j=1, and i-j=0, and MATLAB throws an error.

[Elton Sim]

Hi Ameer, thanks for the response. Do you know if there any alternatives or recommendations to compute the results? This is an assignment in which I am supposed to use that function

[Walter Roberson]

You should probably be ensuring that i < j or j < i

[Elton Sim]

Hi Walter, I can't do that because although it resolves the problem, the results will not be the one I want. There must be cases where i=j for my numerical solution.

[Ameer Hamza]

It is just one case. You will also get zero in the denominator at j=3, and also at i+j-6=0. I don't think there is an easy way to avoid this other than adding all these conditions in for-loop. For example,

for si = 1:n
    for sj = 1:n
        if ~((i==j) | (j==3) | ((i+j-6)==0))
            Mc(si, sj) = subs(sigma, {i,j}, {si,sj});
        end
    end
end