Solving Jacobian matrix for singularities

Hello everyone!

In our assignment we have been tasked with deriving a Jacobian with the help of a provided DH table for a 6DOF robot (UR5). I'm pretty sure I managed to create my jacobian. The next task is specified as follows:

"Analyze the robot for any singularities and if any identify at what joint positions (angles) the singularities occur."

So after a few youtube videos I have concluded that you get the singularities when the matrix for the Jacobian = 0, in other words J[6x6]=0. So my first question to you is if I have understood it correctly.

[Edit] My assumption was incorrect, it's the determinant of the Jacobian that has to equal 0, in other words: det(J[6x6])=0.

If I have understood it correctly, I believe that I can use a solver where I input my my equation as J[6x6]=0 where my variables are TH1, TH2, TH3, TH4, TH5 and TH6. My second question to you would be how I would go about creating a solver for this.

I'm quite new to this so please forgive me if I'm making no sense. For clarification, my Jacobian matrix is entirely symbolic, that means every element is described with some combination of cos(THn), sin(THn) and constants which I have defined as symbolic to keep the matrix somewhat readable.

syms TH1 TH2 TH3 TH4 TH5 TH6
syms d1 d4 d5 d6
syms a2 a3

[d4*cos(TH1) + d6*(cos(TH1)*cos(TH5) + cos(TH2 + TH3 + TH4)*sin(TH1)*sin(TH5)) - a2*cos(TH2)*sin(TH1) - d5*sin(TH1)*sin(TH2 + TH3 + TH4) + a3*sin(TH1)*sin(TH2)*sin(TH3) - a3*cos(TH2)*cos(TH3)*sin(TH1), -cos(TH1)*(a3*sin(TH2 + TH3) - d5*(cos(TH2 + TH3)*cos(TH4) - sin(TH2 + TH3)*sin(TH4)) + a2*sin(TH2) - d6*sin(TH5)*(cos(TH2 + TH3)*sin(TH4) + sin(TH2 + TH3)*cos(TH4))),                           cos(TH1)*(d5*cos(TH2 + TH3 + TH4) - a3*sin(TH2 + TH3) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)),                           cos(TH1)*(d5*cos(TH2 + TH3 + TH4) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)), d6*cos(TH1)*cos(TH2)*cos(TH5)*sin(TH3)*sin(TH4) - d6*cos(TH1)*cos(TH2)*cos(TH3)*cos(TH4)*cos(TH5) - d6*sin(TH1)*sin(TH5) + d6*cos(TH1)*cos(TH3)*cos(TH5)*sin(TH2)*sin(TH4) + d6*cos(TH1)*cos(TH4)*cos(TH5)*sin(TH2)*sin(TH3),                                                            0]
[d4*sin(TH1) + d6*(cos(TH5)*sin(TH1) - cos(TH1)*cos(TH2 + TH3 + TH4)*sin(TH5)) + a2*cos(TH1)*cos(TH2) + d5*cos(TH1)*sin(TH2 + TH3 + TH4) - a3*cos(TH1)*sin(TH2)*sin(TH3) + a3*cos(TH1)*cos(TH2)*cos(TH3), -sin(TH1)*(a3*sin(TH2 + TH3) - d5*(cos(TH2 + TH3)*cos(TH4) - sin(TH2 + TH3)*sin(TH4)) + a2*sin(TH2) - d6*sin(TH5)*(cos(TH2 + TH3)*sin(TH4) + sin(TH2 + TH3)*cos(TH4))),                           sin(TH1)*(d5*cos(TH2 + TH3 + TH4) - a3*sin(TH2 + TH3) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)),                           sin(TH1)*(d5*cos(TH2 + TH3 + TH4) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)), d6*cos(TH1)*sin(TH5) - d6*cos(TH2)*cos(TH3)*cos(TH4)*cos(TH5)*sin(TH1) + d6*cos(TH2)*cos(TH5)*sin(TH1)*sin(TH3)*sin(TH4) + d6*cos(TH3)*cos(TH5)*sin(TH1)*sin(TH2)*sin(TH4) + d6*cos(TH4)*cos(TH5)*sin(TH1)*sin(TH2)*sin(TH3),                                                            0]
[                                                                                                                                                                                                      0,                                      a3*cos(TH2 + TH3) + (d6*sin(TH2 + TH3 + TH4 - TH5))/2 + a2*cos(TH2) + d5*sin(TH2 + TH3 + TH4) - (d6*sin(TH2 + TH3 + TH4 + TH5))/2, a3*cos(TH2 + TH3) + (d6*sin(TH2 + TH3 + TH4 - TH5))/2 + d5*sin(TH2 + TH3 + TH4) - (d6*sin(TH2 + TH3 + TH4 + TH5))/2, (d6*sin(TH2 + TH3 + TH4 - TH5))/2 + d5*sin(TH2 + TH3 + TH4) - (d6*sin(TH2 + TH3 + TH4 + TH5))/2,                                                                                                                                                            -d6*(sin(TH2 + TH3 + TH4 + TH5)/2 + sin(TH2 + TH3 + TH4 - TH5)/2),                                                            0]
[                                                                                                                                                                                                      0,                                                                                                                                                               sin(TH1),                                                                                                            sin(TH1),                                                                                        sin(TH1),                                                                                                                                                                                                cos(TH1)*sin(TH2 + TH3 + TH4),   cos(TH5)*sin(TH1) - cos(TH1)*cos(TH2 + TH3 + TH4)*sin(TH5)]
[                                                                                                                                                                                                      0,                                                                                                                                                              -cos(TH1),                                                                                                           -cos(TH1),                                                                                       -cos(TH1),                                                                                                                                                                                                sin(TH1)*sin(TH2 + TH3 + TH4), - cos(TH1)*cos(TH5) - cos(TH2 + TH3 + TH4)*sin(TH1)*sin(TH5)]
[                                                                                                                                                                                                      1,                                                                                                                                                                      0,                                                                                                                   0,                                                                                               0,                                                                                                                                                                                                        -cos(TH2 + TH3 + TH4),                               -sin(TH5)*sin(TH2 + TH3 + TH4)]
 

If additional information is needed I will provide what I can, thank you for your time!

[Paul]

As to question 1, are you sure the singularities are defined when J = 0? Is it possible you meant that the singularities occur when the determinant of J is zero?

[Daniel Arezoomande]

Hey Paul, thank you for checking in!

Yes so I kept reading up on the subject and you are completely correct, I'll have to correct my post.

So when the det(J[6x6])=0 I get the singularities. Now here is the issue I ran in to. The determinant of the Jacobian gives me a single equation with 6 variables. It is, as far as I know, impossible to calculate my variables unless I have atleast 6 different equations.

Continuing to read up on the subject I found out that I can take determinants of different partial Jacobian so that I get more equations, but so far I have not discovered in what way I can select the partial Jacobian. Would I be able to select which ever elements I want, like 1,1 1,2 1,3 2,1 2,2 2,3 3,1 3,2 3,3?

Answer 1

[Paul]

Becase this is an assignment, here's some code and hints to get started

syms TH1 TH2 TH3 TH4 TH5 TH6
TH = sym('TH',[1 6]); % gather up the unknowns into a single array, makes some things easier
syms d1 d4 d5 d6
syms a2 a3
% define the Jacobian as a single 6 x 6 matrix
J = [
    d4*cos(TH1) + d6*(cos(TH1)*cos(TH5) + cos(TH2 + TH3 + TH4)*sin(TH1)*sin(TH5)) - a2*cos(TH2)*sin(TH1) - d5*sin(TH1)*sin(TH2 + TH3 + TH4) + a3*sin(TH1)*sin(TH2)*sin(TH3) - a3*cos(TH2)*cos(TH3)*sin(TH1), -cos(TH1)*(a3*sin(TH2 + TH3) - d5*(cos(TH2 + TH3)*cos(TH4) - sin(TH2 + TH3)*sin(TH4)) + a2*sin(TH2) - d6*sin(TH5)*(cos(TH2 + TH3)*sin(TH4) + sin(TH2 + TH3)*cos(TH4))),                           cos(TH1)*(d5*cos(TH2 + TH3 + TH4) - a3*sin(TH2 + TH3) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)),                           cos(TH1)*(d5*cos(TH2 + TH3 + TH4) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)), d6*cos(TH1)*cos(TH2)*cos(TH5)*sin(TH3)*sin(TH4) - d6*cos(TH1)*cos(TH2)*cos(TH3)*cos(TH4)*cos(TH5) - d6*sin(TH1)*sin(TH5) + d6*cos(TH1)*cos(TH3)*cos(TH5)*sin(TH2)*sin(TH4) + d6*cos(TH1)*cos(TH4)*cos(TH5)*sin(TH2)*sin(TH3),                                                            0;
    d4*sin(TH1) + d6*(cos(TH5)*sin(TH1) - cos(TH1)*cos(TH2 + TH3 + TH4)*sin(TH5)) + a2*cos(TH1)*cos(TH2) + d5*cos(TH1)*sin(TH2 + TH3 + TH4) - a3*cos(TH1)*sin(TH2)*sin(TH3) + a3*cos(TH1)*cos(TH2)*cos(TH3), -sin(TH1)*(a3*sin(TH2 + TH3) - d5*(cos(TH2 + TH3)*cos(TH4) - sin(TH2 + TH3)*sin(TH4)) + a2*sin(TH2) - d6*sin(TH5)*(cos(TH2 + TH3)*sin(TH4) + sin(TH2 + TH3)*cos(TH4))),                           sin(TH1)*(d5*cos(TH2 + TH3 + TH4) - a3*sin(TH2 + TH3) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)),                           sin(TH1)*(d5*cos(TH2 + TH3 + TH4) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)), d6*cos(TH1)*sin(TH5) - d6*cos(TH2)*cos(TH3)*cos(TH4)*cos(TH5)*sin(TH1) + d6*cos(TH2)*cos(TH5)*sin(TH1)*sin(TH3)*sin(TH4) + d6*cos(TH3)*cos(TH5)*sin(TH1)*sin(TH2)*sin(TH4) + d6*cos(TH4)*cos(TH5)*sin(TH1)*sin(TH2)*sin(TH3),                                                            0;
    0,                                      a3*cos(TH2 + TH3) + (d6*sin(TH2 + TH3 + TH4 - TH5))/2 + a2*cos(TH2) + d5*sin(TH2 + TH3 + TH4) - (d6*sin(TH2 + TH3 + TH4 + TH5))/2, a3*cos(TH2 + TH3) + (d6*sin(TH2 + TH3 + TH4 - TH5))/2 + d5*sin(TH2 + TH3 + TH4) - (d6*sin(TH2 + TH3 + TH4 + TH5))/2, (d6*sin(TH2 + TH3 + TH4 - TH5))/2 + d5*sin(TH2 + TH3 + TH4) - (d6*sin(TH2 + TH3 + TH4 + TH5))/2,                                                                                                                                                            -d6*(sin(TH2 + TH3 + TH4 + TH5)/2 + sin(TH2 + TH3 + TH4 - TH5)/2),                                                            0;
    0,                                                                                                                                                               sin(TH1),                                                                                                            sin(TH1),                                                                                        sin(TH1),                                                                                                                                                                                                cos(TH1)*sin(TH2 + TH3 + TH4),   cos(TH5)*sin(TH1) - cos(TH1)*cos(TH2 + TH3 + TH4)*sin(TH5);
    0,                                                                                                                                                              -cos(TH1),                                                                                                           -cos(TH1),                                                                                       -cos(TH1),                                                                                                                                                                                                sin(TH1)*sin(TH2 + TH3 + TH4), - cos(TH1)*cos(TH5) - cos(TH2 + TH3 + TH4)*sin(TH1)*sin(TH5);
    1,                                                                                                                                                                      0,                                                                                                                   0,                                                                                               0,                                                                                                                                                                                                        -cos(TH2 + TH3 + TH4),                               -sin(TH5)*sin(TH2 + TH3 + TH4)];
theta = sym('theta',[1 6]);
J = subs(J,TH,theta); % use theta to make the output use greek letter symbols
assume(theta,'real'); % the angles are real, so let's assume that to be the case
D(theta) = simplify(det(J),100) % compute the determinant

D(theta1, theta2, theta3, theta4, theta5, theta6) = 

%  a trivial solution for D == 0, is a2 = 0 or a3 = 0, so let's assume those are not true 
assume(a2 ~= 0);
assume(a3 ~= 0);

% now need to solve D(theta) == 0, note: there is one obvious solution
% check out the solve() function, and its option ReturnConditions
% doc solve

[Daniel Arezoomande]

Hey again Paul,

Sorry for bothering you again, but I tried checking out the solve() function with the ReturnConditions, what I gathered from it is that the function will look like this:

X = Dtheta == 0;
THvarlist = sym('THvar',[1 6]);
[THvarlist,parameters,conditions] = solve(X,theta,'ReturnConditions',true)

Where theta is as defined earlier on ,theta = sym('theta', [1 6]), Dtheta is where my determinant of the J function is located.

I'm getting an error message that the number of output variables(1) does not match the input(6). But that's strange because I defined the output variable as a 1 by 6 matrix?

[Paul]

I don't think you can use THvarlist as that output argument the way you'r envisioning. I think you need to specify each of the 6 output solution variables separately, which is a pain.

Referencing back to the code I posted ... it may be easier to return the solution in a single structure:

sol = solve(D==0,theta,'ReturnConditions',true)

The solutions, parameters, and conditions will be fields of sol, which is nice and compact. Try that and feel free to come back and post the results with additional questions if you have any.

[Daniel Arezoomande]

Thank you for being so helpful!

I ran the modified solve function and actually got some "values" for my theta's, which is a marked improvement! I'm a little confused as to what they mean. So from what I gathered, when the det(J) = 0 I'm suppose to have a few different combinations of theta1 theta2 and so on that will fulfill the equation. But the values I'm getting back from the function seem to be variables aswell?

>> sol.theta1
 
ans =
 
 w
z1
z1
z1

>> sol.theta2
 
ans =
                                                                                                                                                                                                                                                        v
2*pi*k - log(((7851*exp(z2*1i) + exp(z2*1i)*exp(z3*1i)*1893i + 8500)*(7851*exp(z3*1i) + 8500*exp(z2*1i)*exp(z3*1i) - 1893i))^(1/2)/(exp((z2*1i)/2)*exp((z3*1i)/2)*8500i + exp((z3*1i)/2)*exp((z2*3i)/2)*7851i - 1893*exp((z2*3i)/2)*exp((z3*3i)/2)))*1i
2*pi*k - log(-((7851*exp(z2*1i) + exp(z2*1i)*exp(z3*1i)*1893i + 8500)*(7851*exp(z3*1i) + 8500*exp(z2*1i)*exp(z3*1i) - 1893i))^(1/2)/(exp((z2*1i)/2)*exp((z3*1i)/2)*8500i + exp((z3*1i)/2)*exp((z2*3i)/2)*7851i - 1893*exp((z2*3i)/2)*exp((z3*3i)/2)))*1i
                                                                                                                                                                                                                                                       x

Am I doing something wrong or are these type of answers what you expect?

[Daniel Arezoomande]

Another thing I'm curious about is

%  a trivial solution for D == 0, is a2 = 0 or a3 = 0, so let's assume those are not true 
assume(a2 ~= 0);
assume(a3 ~= 0);

as my a2 a3 d1 d4 d5 d6 are all constants, but when the values are defined my matrixes are incredibly clutered so I choose to keep them as symbols until I have a function for calculating my thetas. At which point I define my variables as

a2=-0.425;
a3=-0.39255;
d1=0.089159;
d4=0.10915;
d5=0.09465;
d6=0.0823;

[Paul]

It looks like you didn't include the line:

assume(theta,'real')

which was included in the code in my original answer. If you don't include that, Matlab's default assumption is that the thetas can be complex, which makes the solution more .... complex. So make sure to include that assume() statement in your code before the call to solve.

The new values that you are getting back from sol aren't really variables. They are free parameters in the solution that might be subject to additional constraints. This happens when the solutions aren't fixed. For example, consider this simple problem of solving for x1 and x2

syms x1 x2 real
syms y(x1,x2)
y(x1,x2) = x1*sin(x2)

y(x1, x2) = 

The equation y(x1,x2) = 0 should have two solutions: either x1 is 0, in which case the value of x2 can be anything, or x2 is an integer multiple of pi, in which case x1 can be anything. Let's see what solve() returns.

sol = solve(y(x1,x2)==0,[x1 x2],'ReturnConditions',true);
[sol.x1 sol.x2]  % note: sol.x1 and sol.x2 are column vectors

ans = 

sol.parameters

ans = 

sol.conditions

ans = 

You can interpret this solution as follows: sol has two solutions.

In the first solution x1=0 and x2 = z, where z is a free parameter (as shown in sol.parameters). There are no constraints on this solution, as shown by sol.conditions(1) = symtrue.

In the second solution, x1 = z, which is a free parameter, and x2 = k*pi. k is also a free parameter, but as shown in sol.conditions(2), k must be an integer (that's what that funny looking Z means).

So in your problem, solve() may return more than one form of the solution, and each of those could have free parameters, but with different constraints, as shown in this simple example.

I'm not sure what you're asking about in this comment. So I'll just say that if you solve the problem with the symbolic constants, but you don't put the restrictions on a2 and a3 both not being zero, then your symbolic solution will be more complex than it needs to be.

[Daniel Arezoomande]

Wow, not only do you help with coding but you even help me understand the theory behind it, you have gone above and beyond for a random stranger, thank you!

So as usual when I don't quite understand what's happening I try to change some things in my code to see how the output varies. I put my code back in order now

J = simplify([
    Jv
    Jw
    ],100);


assume(theta,'real');

Dtheta = simplify(det(J),100)
assume(a2 ~= 0);
assume(a3 ~= 0);


sol = solve(Dtheta == 0,theta,'ReturnConditions',true)

Your simple example really helped me understand how the solution for the variable can be defined as a parameter and the conditions show me what does parameter values are allowed to be.

>> sol.parameters 
ans =
[k, u, v, w, x, y, z1, z2, z3, z4, z5]

sol.conditions
ans =
        in(k, 'integer') & in(u, 'real') & in(v, 'real') & in(w, 'real') & in(x, 'real') & in(y, 'real')
in(x, 'real') & in(z1, 'real') & in(z3, 'real') & in(z4, 'real') & in(z5, 'real') & in(z2/pi, 'integer')

sol.theta1
ans =
 w
z1

so from my parameters I can see that both w & z1 is a free parameter that has to be an integer that has no complex component.

But I have trouble understanding how I would present my answer. The singularities for my 6DOF robot occur when w and z1 are any real numbers, but would that not mean that I get a singularity for every possible theta1?

[Paul]

Can you execute and post the results of the following commands:

[sol.theta1 sol.theta2 sol.theta3 sol.theta4 sol.theta5 sol.theta6]
sol.parameters
sol.conditions

[Daniel Arezoomande]

Will do!

>> [sol.theta1 sol.theta2 sol.theta3 sol.theta4 sol.theta5 sol.theta6]
sol.parameters
sol.conditions
 
ans =
 
[ w, v,  u,  y, pi*k,  x]
[z1, x, z2, z3,   z4, z5]
 
 
ans =
 
[k, u, v, w, x, y, z1, z2, z3, z4, z5]
 
 
ans =
 
        in(k, 'integer') & in(u, 'real') & in(v, 'real') & in(w, 'real') & in(x, 'real') & in(y, 'real')
in(x, 'real') & in(z1, 'real') & in(z3, 'real') & in(z4, 'real') & in(z5, 'real') & in(z2/pi, 'integer')
 

[Paul]

Excellent.

The first solution says that theta5 must be a multiple of pi and that multiple (k) must be an integer (in(k,'integer)) , and the other thetas are unconstrained real numbers.

The second solution says that theta3 = z2, with the constraint that z2/pi is an integer, and the other thetas are unconstrained. So theta3 = z2 = pi*integer. I don't know why the second solution wasn't returned with theta3 = k*pi like the first solution.

Let's put the whole code together:

syms TH1 TH2 TH3 TH4 TH5 TH6
TH = sym('TH',[1 6]); % gather up the unknowns into a single array, makes some things easier
syms d1 d4 d5 d6
syms a2 a3
% define the Jacobian as a single 6 x 6 matrix
J = [
    d4*cos(TH1) + d6*(cos(TH1)*cos(TH5) + cos(TH2 + TH3 + TH4)*sin(TH1)*sin(TH5)) - a2*cos(TH2)*sin(TH1) - d5*sin(TH1)*sin(TH2 + TH3 + TH4) + a3*sin(TH1)*sin(TH2)*sin(TH3) - a3*cos(TH2)*cos(TH3)*sin(TH1), -cos(TH1)*(a3*sin(TH2 + TH3) - d5*(cos(TH2 + TH3)*cos(TH4) - sin(TH2 + TH3)*sin(TH4)) + a2*sin(TH2) - d6*sin(TH5)*(cos(TH2 + TH3)*sin(TH4) + sin(TH2 + TH3)*cos(TH4))),                           cos(TH1)*(d5*cos(TH2 + TH3 + TH4) - a3*sin(TH2 + TH3) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)),                           cos(TH1)*(d5*cos(TH2 + TH3 + TH4) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)), d6*cos(TH1)*cos(TH2)*cos(TH5)*sin(TH3)*sin(TH4) - d6*cos(TH1)*cos(TH2)*cos(TH3)*cos(TH4)*cos(TH5) - d6*sin(TH1)*sin(TH5) + d6*cos(TH1)*cos(TH3)*cos(TH5)*sin(TH2)*sin(TH4) + d6*cos(TH1)*cos(TH4)*cos(TH5)*sin(TH2)*sin(TH3),                                                            0;
    d4*sin(TH1) + d6*(cos(TH5)*sin(TH1) - cos(TH1)*cos(TH2 + TH3 + TH4)*sin(TH5)) + a2*cos(TH1)*cos(TH2) + d5*cos(TH1)*sin(TH2 + TH3 + TH4) - a3*cos(TH1)*sin(TH2)*sin(TH3) + a3*cos(TH1)*cos(TH2)*cos(TH3), -sin(TH1)*(a3*sin(TH2 + TH3) - d5*(cos(TH2 + TH3)*cos(TH4) - sin(TH2 + TH3)*sin(TH4)) + a2*sin(TH2) - d6*sin(TH5)*(cos(TH2 + TH3)*sin(TH4) + sin(TH2 + TH3)*cos(TH4))),                           sin(TH1)*(d5*cos(TH2 + TH3 + TH4) - a3*sin(TH2 + TH3) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)),                           sin(TH1)*(d5*cos(TH2 + TH3 + TH4) + d6*sin(TH5)*sin(TH2 + TH3 + TH4)), d6*cos(TH1)*sin(TH5) - d6*cos(TH2)*cos(TH3)*cos(TH4)*cos(TH5)*sin(TH1) + d6*cos(TH2)*cos(TH5)*sin(TH1)*sin(TH3)*sin(TH4) + d6*cos(TH3)*cos(TH5)*sin(TH1)*sin(TH2)*sin(TH4) + d6*cos(TH4)*cos(TH5)*sin(TH1)*sin(TH2)*sin(TH3),                                                            0;
    0,                                      a3*cos(TH2 + TH3) + (d6*sin(TH2 + TH3 + TH4 - TH5))/2 + a2*cos(TH2) + d5*sin(TH2 + TH3 + TH4) - (d6*sin(TH2 + TH3 + TH4 + TH5))/2, a3*cos(TH2 + TH3) + (d6*sin(TH2 + TH3 + TH4 - TH5))/2 + d5*sin(TH2 + TH3 + TH4) - (d6*sin(TH2 + TH3 + TH4 + TH5))/2, (d6*sin(TH2 + TH3 + TH4 - TH5))/2 + d5*sin(TH2 + TH3 + TH4) - (d6*sin(TH2 + TH3 + TH4 + TH5))/2,                                                                                                                                                            -d6*(sin(TH2 + TH3 + TH4 + TH5)/2 + sin(TH2 + TH3 + TH4 - TH5)/2),                                                            0;
    0,                                                                                                                                                               sin(TH1),                                                                                                            sin(TH1),                                                                                        sin(TH1),                                                                                                                                                                                                cos(TH1)*sin(TH2 + TH3 + TH4),   cos(TH5)*sin(TH1) - cos(TH1)*cos(TH2 + TH3 + TH4)*sin(TH5);
    0,                                                                                                                                                              -cos(TH1),                                                                                                           -cos(TH1),                                                                                       -cos(TH1),                                                                                                                                                                                                sin(TH1)*sin(TH2 + TH3 + TH4), - cos(TH1)*cos(TH5) - cos(TH2 + TH3 + TH4)*sin(TH1)*sin(TH5);
    1,                                                                                                                                                                      0,                                                                                                                   0,                                                                                               0,                                                                                                                                                                                                        -cos(TH2 + TH3 + TH4),                               -sin(TH5)*sin(TH2 + TH3 + TH4)];
theta = sym('theta',[1 6]);
J = subs(J,TH,theta); % use theta to make the output use greek letter symbols
assume(theta,'real'); % the angles are real, so let's assume that to be the case
D(theta) = simplify(det(J),100) % compute the determinant

D(theta1, theta2, theta3, theta4, theta5, theta6) = 

%  a trivial solution for D == 0, is a2 = 0 or a3 = 0, so let's assume those are not true 
assume(a2 ~= 0);
assume(a3 ~= 0);
sol = solve(D==0,theta,'ReturnConditions',true);
[sol.theta1 sol.theta2 sol.theta3 sol.theta4 sol.theta5 sol.theta6]

ans = 

sol.parameters

ans = 

sol.conditions

ans = 

The first solution is self evident because of the sin(theta5) out front in the expression for D(theta)

The second solution is not as obvious, be we can verify it

syms k integer
simplify(D(TH1,TH2,k*sym(pi),TH4,TH5,TH6))

ans = 

0

It appears that the solution is any set of thetas (within their allowable physical range) where theta3 = k*pi or theta5 = k*pi, where k is an integer.

The next very important question you should ask yourself is: does this solution make sense based on the physical configuration of your application? If not, it probably indicates an error in forming the Jacobian in the first place.

[Daniel Arezoomande]

It all makes sense now. So if we assume my Jacobian was derived correctly (which I'll be quadruple checking), my conclusions is that the singularities for my robot occur when theta3 or theta5 approaches the angle k*pi regardless of what the other angle values are.

From what I gather the reason we want to know this is so we can add an exclusion for our robot that if theta3 or theta5 approaches k*pi, the robot should look for a different configuration to approach the desired location. Because even though the angles arent exactly k*pi, as we get closer to the singularities, the calculations for linear/angular velocities approach infinity which isn't good either. This I suppose is no longer relevant to the question I posed at the beginning, it's just my own curiosity and willingess to understand robotics.

Regardless, thank you Paul, and thank you for sticking with me. If I find that my Jacobian is incorrect I'll be able to correct it, derive my singularities and understand it all without bothering anyone haha. Lets see if I can manage to solve my last problem in the assignment with the revelations I got from here!

[Paul]

HTH. Good luck.

[Daniel Arezoomande]

As I was explaining my answer I thought back to the sin(theta5) outside of the expression which has an obvious answer. Why does the code not give us a second column of answers where theta5 = k*0?

Because my answer will have to include that singularities occur for when theta3=pi*k, theta5=pi*k,theta5=0 right?

[Daniel Arezoomande]

Nevermind, I guess it's not needed because the case pi*k for when k = 0 gives us 0. Which means theta5=pi*k includes our case when theta5=0.