Strange simulink error. Professor & TA don't know. Any Ideas?

Question

I am currently working on developing a simulation for a restricted UAV, and I’m having some issues with Simulink while doing this. I have my mathematical equations which have been approved by my professor, but while I’m trying to simulate the system I am running into strange issues. Moreover, neither my professor nor my TA can figure out what is causing this issue, despite multiple private sessions to troubleshoot. If anyone can shed some light on the issue it would be highly appreciated.

Firstly, the equation I’ve modeled on Simulink is this:

My Simulink block diagram is above. Input 1 Is Fb+Ff which, in this scenario = 6.8. Input 2 is theta which = 0.

I have annotated the values that should be and are at each of the nodes before the summing junction. In this section of the block diagram the theoretical value and the simulation agree.

The problem comes in the more complex half of the diagram. The graph I am validating should look like this.

There is an initial rise and then oscillations about a non-zero number. Obviously, my graph will be ideal so oscillations will not persist like they do in the graph above. My final graph is this…

Where the y-axis is elevation angle in degrees and the x-axis is in seconds. Obviously these two graphs are in no way similar. In an attempt to analyze where I could have made an error in my block diagram, I analyzed the behavior at each node.

Immediately after the cos block to the lower right I get oscillating behavior

Immediately after the summing junction and gain (1/0.91) I get:

The above bears a resemblance to the graph I’m validating against.

After the first integrator, however, I get this.

Adding the two graphs together, should I not be getting a graph that rises and then oscillates about a non-zero number?

Have I made a mistake in my block diagram somewhere?

The only value I’m allowed to change is the gain marked De.

At this point I am starting to hit a wall and really need some insight on this. Surely this cant be something wrong with the Simulink source code.

Nothing I do gets me anywhere near the validation graph. Can someone please point me in the right direction or let me know if you see any errors?

Any insight is highly appreciated, thank you for reading.

Should you want to view the Simulink file for further analysis, it is available for download here: https://ufile.io/08aqv However please note that the attached file is my entire Simulink project. The circuit reference in the question above is contained in the subsystem "elevation".

Answer 1

Your results are correct. I have solved the differential equations in Octave using the ode solver and I get the same results:

function dX = elevation(t,X)
  g = 9.81;
  m_cw = 1.87;
  l_w = 0.47;
  m_uav = 1.17;
  l_a = 0.66;
  J_e = 0.91;
  D_e = 0.04;
  F_in = 0.68;
  theta = 0;

  dX(1) = X(2);
  dX(2) = (F_in * l_a * cos(theta) + (m_cw*l_w - m_uav*l_a)*g*cos(X(1)) - D_e * X(2)) / J_e;

Which gets saved as elevation.m and solved by:

tspan = [0 60];
X0=[0,0];
opts = odeset('MaxStep',0.1,'InitialStep',0.1);
[t,X] = ode45(@elevation,tspan,X0,opts);
plot (t,X(:,1))
grid
xlabel('time')
ylabel('elevation [rad]')
title ('elevation vs. time')

and which gives the following results:

So the only possible explanations are:

1. The equations you use are incorrect (unlikely because you say they have been checked)
2. The initial conditions of the differential equations are incorrect or do not match those used to get the validated results (assumed to be zero for both elevation and elevation rate at t=0). For example, it looks like the initial condition for the elevation in the validated graph is about -30 deg or so.
3. The numerical values for the various parameters are incorrect or do not match those used to get the validated results
4. The inputs used (Fb+Ff and theta) do not match those used to get the validated results.

To me, #4 seems like the most plausible explanation. You have assumed both inputs to be constant. Is that really the case? Are they not a function of time?

I have tried changing the value of D_e as you suggest and you do get some quite different results. For example, this is with D_e = 0.2 instead of 0.04: