Changing delay length to solve algebraic loops

Edit: Uh oh. I found out that it worked because my signals were one seconds long and was exactly starting at for example sec 1 and ending at sec 3. After changing to non-integer time values, that method also does not works. But still wonder why.

I am using an algebraic loop for a mass-change model. Basically it subtracts the mass of fuel used from the initial mass of vessel. The initial mass is 25 kg and the subtracted mass is variable with time and throttle. At first I used an initial condition block for the process, but eventually it just crashes with an algebraic loop error.

So i tried to use an delay block for the job. It works but the problem is, the values are totally wrong when i use it with a delay length of 1:

I included my mdot mass flow rates to better explain my range of subtraction. Certainly, the output values are totally wrong as i am applying them for just a few seconds. So, i tried to increase the delay length to 1000, and the results are almost identical to what it is supposed to be:

There is some peaks at the end, but i think they are not the main problem here. I wonder why when the delay increases, the outputs get more realistic?

Also, is this a correct way to solve an algebraic loop, at least for my mass change model?

[Paul]

Why is there a loop at all? Isn't the equation just

m(t) = M0 - massused(t)

[astatum]

@Paul It needs to be real-time as the controller response and fuel use are different every run. By the way, I really am not sure how should I do this. Isn't what I done and your equation the same?

Answer 1

[Paul]

Refering to the block diagram at the top of the question ....

Would it be correct to replace the IC block with a Constant block with the "Constant value" parameter set to 25? There is no input ot the Constant block, hence no loop.

[astatum]

Yeah but since I want to dynamically change the mass, its not possible with no loop. I am trying to extract the used fuel mass from the current total mass every step, not the initial mass.

[Paul]

Whether or not a loop is needed depends on the equation for the mass that we're trying to model. Presumably, mass changes as function of time continously, In which case two options are:

1. M(t) = M0 - massused(t), where massused(t) is the total mass expended from t = 0 and can be computed with a closed form expression

2. Mdot(t) = f(t), M(0) = M0, i.e., a differential equation that has to be integrated numerically becuase even though f(t) can be computed, we don't have a close from expression for its anti-derivative.

Or, I suppose you could discretize and only update the mass at discrete times, so something like

3. M(dt*(k + 1)) = M(dt*k) + dt*f(dt*k)

Equation 3 sounds like what you might be trying to do, except that dt*f(dt*k) is the increment to used fuel mass, not the total used fuel mass. But it's hard to say w/o knowing the definitions of t_peri, t_cent, and t_RCS

If you write out the mass equation mathematically, we can probably figure out how to implement it correctly in Simulink.

[astatum]

Alright, my approach to the problem was like this:

I know the mass flow rates of my three different thrusters: central, peripheral and RCS thrusters. Let's say the central thruster got mass flow rate of 0.15 kg/s, the peripheral got 0.3 kg/s, RCS got 0.05 kg/s.

My logic for the solution is to simply get an active-time of the thruster signal (t_peri, t_cent, t_rcs) from the controller block and just multiply it with mass flow rates. After summing them all, extracting it from the current mass (initial mass for the first operation, of course) would give me the live mass. I'm saying live mass because I am going to feed it into the equations of motion block during ongoing simulation.

I will probably screw up with the equation but my best try is here:

for any of the thrusters: m_live = m_live - (mdot*time_on)

So yeah, got an algebraic loop here too. Guess I need to get better at this.

I think that the Equation 3 looks really promising, since it is not a necessity to use the total used fuel mass to change the final value. Since it just works a little different to achieve the same thing, it can solve my problem. I'll try it and reply back, many thanks!

[Paul]

Because m_live is going to the equations of motion, it is crucial that it have no discontinuities. Otherwise, you'll have discontinuous jumps in the vehicle acceleration that are not physically possible.

So the next question is: how are the "time_on" variables computed for each thruster?

I'm going to take a guess that the simulation doesn't really need to compute "time_on" for each thruster. Rather, there's some logic (or a controller) that determines at each instant in time whether a thruster is on or off. In this case, at each instant in time, for each thruster we have a variable "thruster_is_on" which takes on the value 1 or 0. So at any instant in time.the mdot for the central thruster is

mdot_central = 0..15*central_thruster_is_on,

and then the same for the peripheral and the RCS. The total m_dot is

mdot = mdot_central + mdot_peri + mdot_RCS

Now, feed mdot into an Integrator block, and set the inital condition parameter in the integrator block to the initial mass (25 kg).

The output of the integrator is m_live.

Make sure to simulate with a variable step solver and have the zero crossing detection turned on (I believe it is by default) so the simulaton catches the jumps in the thruster on/off logic.

[astatum]

Since I am using a discrete-time solver for simulink due to my project conditions, I just tried it with discrete time integrator. And it works :))

Here is a picture with real actuation time inputs fed into the system:

I used an first-order hold block to make sure that signal is pseudo-continuous, for the sake of smooth equations of motion calculation.

Many thanks!