Thoughts on overhauling the MATLAB S-Function?

[noah-curran]

Hi all, thank you for your hard work on this project. After a long search for a good tool for data collection for my current research paper, we landed on JSBSim after discovering the MATLAB S-Function that @sthelia put together (thanks a ton!).

I've noticed it could use some major work with cleaning up and trimming some unnecessary parts / bugs in it. For instance, the altitude outputs to meters instead of feet (I think this is a bug), and adding documentation for MacOS compiling (just worked this out this past weekend).

I'm wondering if there is any interest in opening up a discussion if I were to overhaul some parts of it? This project seems promising for my own research, so if I were to make some changes while conducting my research I'd love to discuss contributing what I can.

Thanks,
Noah

[seanmcleod]

Out of interest what is your research topic?

I'd say continue the discussion here in terms of your proposed improvements. I did discuss 1 or 2 proposed changes to the S-Function with @agodemar who wrote the original S-Function many years ago, but I never got round to making them.

[noah-curran]

Our research requires studying the effects of anomalous/adversarial control inputs from both autonomous and manual systems. The S-Function gives us the ability to create a control loop that we can inject these anomalies/attacks into.

I'll draft a few of the proposed improvements when I get into the office tomorrow. It's possible I misunderstand parts of the code and what I'm thinking is already possible, unnecessary, or outside of the scope of this S-Function.

[agodemar]

Hi Noah, if you have spotted bugs or possible improvements in the S-Function code, we can definitely work on it.
In my opinion the S-Function should be cleaned up in terms of code style. The design can be also improved, allowing a more general approach in terms of inputs. Also, the possibility to inject winds and turbulence should be much appreciated by users.

[noah-curran]

@agodemar I was actually thinking of the wind injection stuff, too! We may need it for our research, so I was going to look into it.

[seanmcleod]

@noah-curran here is a link to #498 (comment) some of my comments about modifying the current s-function.

In terms of your setup, I presume the control loop, say for example one implementing a g-command FBW law, is going to be implemented in Simulink? As opposed to being implemented in the JSBSim aircraft model? Given that both options are available when using the JSBSim s-function.

If you're modelling actuators, again will they be modelled in Simulink which then simply outputs a control surface deflection angle to JSBSim, or will you use the actuator model available in JSBSim?

Which leads me to my final question. Are the anomalous/adversarial control inputs going to be input at the level of the command law, e.g. if 4g was commanded, will the anomalous control input be modelled as changing that 4g to some delta, or will the anomalous inputs be implemented by adding some delta to the final control surface deflection?

[noah-curran]

Thanks for that link. The discussion was rather helpful and confirms some of my findings from nose diving into the code in the last week.

As for your questions, the short answer is that I believe your assumptions of what I'm trying to do are mostly correct. I think you have the general idea of what I'm trying to achieve, but I'll describe more below.

Here is the longer answer:

1. Yes, we want to implement some control loop in Simulink, similar to how you describe.
2. I think for the actuators it would be better for those to be modeled in JSBSim. Though, I don't know if it is possible to do this given that we want to take control of implement the control loop in Simulink. It currently seems that there is an all or nothing approach happening with respect to loading in scripts (i.e., the option for using the script in the S-Function's block parameter is enabled or not), and I'm not sure if it's possible to pick and choose in this way.
3. And kinda, the anomalous control input would essentially be modeled based on anomalous observations of the state. So suppose some control input is supposed to engage when event E occurs in the state of the airplane. Now, suppose event E doesn't occur in the physical world, but the sensors think it does occur. Thus, the control input will engage following the anomaly that the system believes E is happening, and it follows that an anomalous input to the control system occurs.

Pardon me if I use any terminology that is inconsistent with aircraft/control theory. I come from a computer science background and have been self learning much of this stuff in the past half year, so I'm still pretty rough on the edges. Below is a rough figure to give you an idea of a visual in case there is any confusion in my communication. I tried to generalize it as much as possible, but obviously the generalization may not hold in all cases.

[agodemar]

@noah-curran @seanmcleod
The following is more or less a reproduction of the above scheme with the new tools available in GitHub's rendering engine, that is, Mermaid and MathJax.

  graph TD;
      A[Weather Turbulence] --> C((Control<br/>Decision-Maker));
      B[Pilot Override Input] -->C;
      C -- "Control input" ----> D(JSBSim);
      D -- "Outputs of<br/>aircraft state,<br/>controls,<br>pilot, etc." ----> E>Inject Anomalies];
      E -- "Feedback" ----> C;

Loading

Aircraft state vector:
$$ \boldsymbol{x} = \left[\; x_1, \; x_2, \; \ldots \;, \; x_{N_x} \;\right]^\mathrm{T} $$

Control input vector:
$$ \boldsymbol{u} = \left[\; u_1, \; u_2, \; \ldots \;, \; u_{N_u} \;\right]^\mathrm{T} $$

Pilot, auxiliary variables, etc.:
$$ \boldsymbol{p} = \left[\; p_1, \; p_2, \; \ldots \;, \; p_{N_p} \;\right]^\mathrm{T} $$

Outputs:
$$ \boldsymbol{y} = \big[ \; \boldsymbol{x}; \; \boldsymbol{u}; \; \boldsymbol{p} \; \big] $$

Feedback to Control Decision-Maker block:
$$ \boldsymbol{f} \; = \; \text{outputs} \; + \; \text{anomalies} \; = \; \boldsymbol{y} + \Delta\boldsymbol{y} $$

[noah-curran]

Interesting tool! Thanks for sharing.

[noah-curran]

I'm going to be doing development of the MATLAB S-Function on my fork of JSBSim:
https://github.com/noah-curran/jsbsim

A few things I've noticed and made some changes to:

1. Altitude was in meters instead of feet (labeling indicates it should be in feet)
2. There were a few unused/redundant functions. I either cleaned these up or removed them.
3. Some formatting/spacing.
4. Constants for the sizes of the input/output parameters of the S-Function.
5. Added an additional engine to control in the case of aircraft with two engines.

In making these changes, I'm starting to appreciate @agodemar's suggestion to generalize the parameters. Besides for the problem of the number of engines capable of being controlled, I think there are several other control inputs that could become necessary in the future. Some other state output could also become useful. While it might be a few days of labor at least, I think generalizing it will just make it less painful. I'm probably going to dig into this tomorrow and see where it goes.

[bcoconni]

@noah-curran, FWIW do not hesitate to submit pull requests of small updates as you are progressing in your work. These would be easier to review rather than a big chunk of modified code.

[noah-curran]

@noah-curran, FWIW do not hesitate to submit pull requests of small updates as you are progressing in your work. These would be easier to review rather than a big chunk of modified code.

Sounds good @bcoconni. That's probably the better way to go. I just finished up my first large round of changes late last night, so I'll submit a PR today for y'all and we can iterate from there. Going to do some different development in the coming days and will be starting an internship next week, so apologies if I don't get back around to being able to commit much time to this anytime soon.

I'm also having an undergrad do some QA with my fork, so we'll have some additional eyes on it as well.

[seanmcleod]

In terms of generalizing this is what I was suggesting.

What I was wondering while glancing at the source code was whether you couldn't make the whole set of states and output ports dynamic in terms of reading a set of property names from a config file and use that to set up the output ports during initialization etc. and also use it during each time step to retrieve the property values from JSBSim.

That way any user who wanted to use the JSBSim S-Function with their own unique set of states could simply edit the config file and therefore wouldn't need to understand any C++ coding and wouldn't need to have to recompile the S-Function DLL.

So some sort of ini file along these lines.

[OutputPort 0]
velocities/u-fps
velocities/v-fps
velocities/w-fps
velocities/p-rad_sec
velocities/q-rad_sec
velocities/r-rad_sec
position/h-sl-ft
position/long-gc-rad
position/lat-gc-rad
attitude/phi-rad
attitude/theta-rad
attitude/psi-rad

[OutputPort 1]
...
...

So if I'm working with a piston engined aircraft then I would include properties for mixture, but exclude them if I'm working on a turbine engine model. If I'm interested in fuel flow for example then I include the relevant property to have it available for my Simulink model otherwise I don't bother including it.

If I'm working with a multi-engine model then I include properties for multiple engines, e.g. rpm[0] rpm[1] etc.

As opposed to have some massive union of all/most properties that clutter up the Simulink port blocks most of which I don't want, or if the subset we currently expose via the hard-coded set in the C++ source doesn't include what I need then I have to go and edit C++ code and recompile etc.

[noah-curran]

Thanks for the write-up @seanmcleod. This is something I was planning to work on because the hardcoding/cluttered properties was starting to get pretty annoying. It helps to have a spec for what you guys want, so this will make it a lot easier for me to put together, I think.

One question, though. Currently, the S-Function calls functions such as those in FGFCS.h (e.g., SetDeCmd(double)). Do you think it would be fine to go the route of just directly setting the value through SetPropertyValue(char*, double)? It would surely make things a lot easier, but if there is some underlying code that makes this the wrong way to go, it would be helpful to know before I start.

[bcoconni]

Regarding the .ini file format, I would suggest that we discuss about using such a format prior to making too much coding. Indeed we should try as much as possible to avoid growing the number of input formats that are used by JSBSim, including its bindings.

I know that XML is losing some speed (to say the least) nowadays while it was all the rush when the development of JSBSim has started in the late 90's, early 2000's. Some other formats (YAML, JSON, etc.) have much more visibility and they might actually be better choices than XML for the particular application at discussion. However, I would like the decision being made collectively before rushing to development.

Just my 2 cents.

[seanmcleod]

I don't have an overall spec. These were just some of my observations when I reviewed the code when it was submitted. I don't actually use Matlab/Simulink for my work, I use JSBSim pretty much as is with some minor tweaks of my own in the VSS (Variable Stability System) simulator that I developed to help with test pilot training.

Do you think it would be fine to go the route of just directly setting the value through SetPropertyValue(char*, double)

Don't see why not.

PropertyManager->Tie("fcs/elevator-cmd-norm", this, &FGFCS::GetDeCmd, &FGFCS::SetDeCmd);

One thing though I would do is cache the property lookup for performance reasons, especially if users want to run thousands of simulations faster than real-time, see @bcoconni's comment - #530 (comment)

FGPropertyManager* pm = fdm->GetPropertyManager();
FGPropertyNode* node = pm->GetNode("/my/property/name");

while (fdm->Run()) {
  x = node->getDoubleValue();
}

So cache the property nodes for each property specified in the Simulink ports.

[agodemar]

Regarding the file that configures the ports, we can adopt the XML format, along the lines:

<?xml version="1.0"?>
<s_function_config>
  <outputs>
    <output name="Output 1">
      <property> velocities/vc-kts </property>
      <property> accelerations/a-pilot-x-ft_sec2 </property>
      <property> accelerations/a-pilot-y-ft_sec2 </property>
      <property> accelerations/a-pilot-z-ft_sec2 </property>
      <property> attitude/phi-rad </property>
      <property> attitude/theta-rad </property>
      <property> attitude/psi-rad </property>
    </output>
    <output name="Output 2">
      <property> position/h-agl-ft </property>
      <property> position/lat-gc-deg </property>
      <property> position/long-gc-deg </property>
    </output>
  </outputs>
</s_function_config>

[noah-curran]

@agodemar @seanmcleod Thanks for this tips. I think this shouldn't be too difficult of a task.

[noah-curran]

@agodemar @seanmcleod Great news, I have a preliminary version of the configurable port input/output feature working pretty well on my fork:
https://github.com/noah-curran/jsbsim/tree/master/matlab

One thing to note is that I also added a <input> property (however, only one input port is allowed) in addition to the <outputs> property. Here is the test xml file I used:

<?xml version="1.0"?>
<s_function_config>
  <input>
    <property> fcs/elevator-cmd-norm </property>
  </input>
  <outputs>
    <output name="Output 1">
      <property> velocities/vc-kts </property>
      <property> accelerations/a-pilot-x-ft_sec2 </property>
      <property> accelerations/a-pilot-y-ft_sec2 </property>
      <property> accelerations/a-pilot-z-ft_sec2 </property>
      <property> attitude/phi-rad </property>
      <property> attitude/theta-rad </property>
      <property> attitude/psi-rad </property>
    </output>
    <output name="Output 2">
      <property>position/h-agl-ft</property>
      <property>position/lat-gc-deg</property>
      <property>position/long-gc-deg</property>
    </output>
  </outputs>
</s_function_config>

Another thing to note is that the ports, while named in the xml file, cannot be dynamically named in the S-Function, as it is unsupported by MATLAB at the moment. Who knows if they will ever add this feature, but if they do I made a comment in the code of where it would be best to add it.

A final thing to note is that because the I/O ports to the S-Function are dynamic, I believe a single initial run is required to get the ports to pop up in the correct configuration. I haven't had enough time to play around with this to confirm this though.

I'm sure there is a bug or two, but I'm going to tackle a stress test tomorrow to double check and do some clean-up.

The next task I'm going to go after is allowing weather input to the S-Function.

[bcoconni]

@noah-curran nice to see that the XML format has been selected 😄. Ignore my comment #627 (reply in thread) above.

[noah-curran]

Ignore my comment #627 (reply in thread) above.

Sounds good. Though, I think even given the other options, XML is still a good choice for this particular problem since we can have attributes. While the attributes aren't currently used, I think they may allow us to have more flexibility with adding future features. (For instance, if MATLAB ever gives us the ability to name ports, we can simply use the "name" attribute to optionally set names for the ports.)

[agodemar]

Great news @noah-curran
Could we consider to have an input signal with a dynamic dimension? That is, it could be a scalar signal (a 1D column of values per each discrete time value) or a multi-column signal (1 column per each input)? In principle, this is possible in standard Simulink workflows.

Regarding the initialization of the simulink block, I would try the callback mechanism. That is, you can use callbacks for various events in the lifecycle of a Simulink project. There is one event that is triggered when the project is loaded up, that we can use to initialize the S-Function block.

[noah-curran]

Are you thinking of something like this:

<?xml version="1.0"?>
<s_function_config>
  <input type="vector|scalar">
    <property> prop1 </property>
    <property> prop2 </property>
    <property> prop3 </property>
  </input>
  <outputs>
   ...
  </outputs>
</s_function_config>

Where the default value of attribute type of the input element could be "scalar" just because that is probably the easier input style for users to use out of the box.

There is one event that is triggered when the project is loaded up, that we can use to initialize the S-Function block.

This could be a good solution for now. One problem I think foresee is that if the Simulink project isn't loaded up and one opts to use sim(...) in a MATLAB script, then the probject is never given a chance to load up. I doubt this is fixable, but I also doubt one would usually run a simulation that they haven't confirmed is working by first opening the Simulink file.

Another problem I could foresee is that if someone inserts this block, sets the parameters, and applies them. Maybe there is a callback for this event that I could make use of as well, so this might be fixable.

[agodemar]

One problem I think foresee is that if the Simulink project isn't loaded up and one opts to use sim(...) in a MATLAB script, then the probject is never given a chance to load up. I doubt this is fixable, but I also doubt one would usually run a simulation that they haven't confirmed is working by first opening the Simulink file.

@noah-curran there are other events that might trigger a callback. These can be captured at different levels: at model level, at block level, or at port level. Here is the online documentation.
For instance, at the S-Function block level we can use the standard events/callbacks named InitFcn or StartFcn.
See this documentation page.

[noah-curran]

@agodemar BTW- the callbacks you linked there aren't the one's the C MEX S-Function uses. You can find the choice of callback functions available to us on this documentation page (be sure to look under the C MEX header). I chose mdlCheckParameters and mdlProcessParameters.

mdlCheckParameters is first called when parameters are changed. This will do error checking to validate that the input looks OK.

mdlProcessParameters is called immediately after in order to process the changes. Note that this is also called by MATLAB when the input is incorrect in order to refresh it to the last correct input. I think they do this in case a developer does not add error handling to this function for when the validation returns an error. I notice it is called twice when the validation is erroneous, so this means it is first called to see if the invalid input will still work and then calls it again to refresh to the last valid input. Because invalid input can be breaking to the processing code, I add error checking anyway. This isn't really important information, but just adding explanation for the MATLAB callback lifecycle.

These functions are also called in mdlInitializeSizes. The MATLAB documentation for mdlCheckParameters recommended this way of organizing the code. Without this call, the S-Function will break on the first run.

I looked through the other callbacks available to us, and these should be the only ones of any use for this task of updating the input/output ports of the JSBSim S-Function block. I've tested it, and each of the problems mentioned before are now resolved. I've committed this code to my fork.

[seanmcleod]

a 1D column of values per each discrete time value

So I think what @agodemar is requesting is an array of values for the input to be fed in on each time step, e.g. say for fcs/elevator-cmd he can supply a pre-set of inputs for this input variable.

I guess for example he could pre-generate a doublet input.

[agodemar]

So I think what @agodemar is requesting is an array of values for the input to be fed in on each time step, e.g. say for fcs/elevator-cmd he can supply a pre-set of inputs for this input variable.
I guess for example he could pre-generate a doublet input.

@seanmcleod Correct. Or, the values could be generated by a controller designed in Simulink.
When the array becomes a multi-column array, each column would represent a separate channel. For instance, elevator command, rudder, speed breaker, etc.

[noah-curran]

@agodemar Hmm... I see what you're getting at, but why is it the S-Function's job to process the input in this way? Source blocks are normally given the task of taking a set of data and feeding it to the Simulink pipeline one data point at a time. For example, the from workspace block.

Is there an example of another block that takes several time steps of data all at once and then processes it like this? Is there a specific feature you're hoping to give the S-Function the capability of having by having more data loaded in all at once? Sorry for all of the questions here, I'm just not understanding the motivation for having this feature. If something else can do this job.

[noah-curran]

@seanmcleod Honestly, I don't see why that would be useful here. Why not just generate your values and feed them to the S-Function one at time? This is what every data source block does in Simulink. I have never seen values stack up in this way besides for in Laplace functions that have a time delay on input. Since in this case there is no need for this functionality for JSBSim's input, I don't see the usefulness here.

[seanmcleod]

Honestly, I don't see why that would be useful here.

@agodemar will have to motivate given his usage of Simulink.

[seanmcleod]

@noah-curran ready to submit a pull request yet?

[noah-curran]

I'll go ahead and submit one. However, when we tried to compile on Windows machines we noticed a few issues and designed a workaround that that worked for us. I'll outline it in the PR.

[noah-curran]

See here: #651