Skip to content

Commit

Permalink
- adding thresholds and gain scheduling inded to gainsched model
Browse files Browse the repository at this point in the history
  • Loading branch information
Steinar Elgsæter committed Nov 1, 2023
1 parent 83b8939 commit b2e6e00
Show file tree
Hide file tree
Showing 4 changed files with 74 additions and 20 deletions.
48 changes: 41 additions & 7 deletions Dynamic/SimulatableModels/GainSchedModel.cs
View file
Original file line number Diff line number Diff line change
Expand Up @@ -241,6 +241,7 @@ public void SetModelParameters(GainSchedParameters parameters)
else*/
{
double x_otherInputs = modelParameters.Bias;
double gainSched = givenInputs[modelParameters.GainSchedParameterIndex];
//nb! input may include a disturbance!
if (givenInputs != null)
{
Expand All @@ -250,7 +251,7 @@ public void SetModelParameters(GainSchedParameters parameters)
continue;
if (i < GetModelInputIDs().Length)//model inputs
{
x_otherInputs += CalculateLinearProcessGainTerm(i, givenInputs[i]);
x_otherInputs += CalculateLinearProcessGainTerm(i, givenInputs[i], gainSched);
}
else // additive inputs
{
Expand Down Expand Up @@ -279,19 +280,50 @@ public void SetModelParameters(GainSchedParameters parameters)
/// </summary>
/// <param name="inputIndex">the index of the input</param>
/// <param name="u">the value of the input</param>
/// <param name="u_GainSched">the value of scheduling input</param>
/// <returns>contribution to the output y, excluding bias and curvature contributions</returns>
private double CalculateLinearProcessGainTerm(int inputIndex, double u)
private double CalculateLinearProcessGainTerm(int inputIndex, double u, double u_GainSched)
{
double processGainTerm = 0;
int gainSchedModelIdx = 0;
for (int idx = 0; idx < modelParameters.LinearGainThresholds.Length; idx++)
{
if (idx == 0)
{
if (u_GainSched < modelParameters.LinearGainThresholds[idx])
{
gainSchedModelIdx = idx;
break;// jump out of for-loop
}
else
{
gainSchedModelIdx = idx+1;
}
}
else if (idx == modelParameters.LinearGainThresholds.Length - 1)
{
if (u_GainSched > modelParameters.LinearGainThresholds[idx])
{
gainSchedModelIdx = idx;
}
}
else
{
if (u_GainSched > modelParameters.LinearGainThresholds[idx-1] &&
u_GainSched < modelParameters.LinearGainThresholds[idx])
{
gainSchedModelIdx = idx;
}
}
}

if (modelParameters.U0 != null)
{
//TODO: use correct linear gain
processGainTerm += modelParameters.LinearGains.First()[inputIndex] * (u- modelParameters.U0[inputIndex]);
processGainTerm += modelParameters.LinearGains.ElementAt(gainSchedModelIdx)[inputIndex] * (u- modelParameters.U0[inputIndex]);
}
else
{
//TODO: use correct linear gain
processGainTerm += modelParameters.LinearGains.First()[inputIndex] * u;
processGainTerm += modelParameters.LinearGains.ElementAt(gainSchedModelIdx)[inputIndex] * u;
}
return processGainTerm;
}
Expand All @@ -309,6 +341,8 @@ private double CalculateStaticStateWithoutAdditive(double[] inputs, double badVa
double x_static = modelParameters.Bias;

// inputs U may include a disturbance as the last entry
double gainSched = inputs[modelParameters.GainSchedParameterIndex];

for (int curInput = 0; curInput < Math.Min(inputs.Length, GetLengthOfInputVector()); curInput++)
{
if (curInput + 1 <= modelParameters.GetNumInputs())
Expand All @@ -326,7 +360,7 @@ private double CalculateStaticStateWithoutAdditive(double[] inputs, double badVa
}
lastGoodValuesOfInputs[curInput] = inputs[curInput];
}
x_static += CalculateLinearProcessGainTerm(curInput, curUvalue);
x_static += CalculateLinearProcessGainTerm(curInput, curUvalue, gainSched);
}
}
return x_static;
Expand Down
19 changes: 17 additions & 2 deletions Dynamic/SimulatableModels/GainSchedParameters.cs
View file
Original file line number Diff line number Diff line change
Expand Up @@ -27,20 +27,30 @@ public class GainSchedParameters : ModelParametersBaseClass

public FittingSpecs FittingSpecs = new FittingSpecs();


/// <summary>
/// A time constant in seconds, the time a 1. order linear system requires to do 63% of a step response.
/// Set to zero to turn off time constant in model.
/// </summary>
public double[] TimeConstant_s { get; set; } =null;



/// <summary>
/// The index of the parameter among the ModelInputIdx
/// </summary>
public int GainSchedParameterIndex = 0;


/// <summary>
/// Thresholds for when to use which gains
/// </summary>
public double[] TimeConstantThresholds { get; set; } = null;

/// <summary>
/// The uncertinty of the time constant estimate
/// </summary>
public double[] TimeConstantUnc_s { get; set; } = null;


/// <summary>
/// The time delay in seconds.This number needs to be a multiple of the sampling rate.
/// Set to zero to turn of time delay in model.
Expand All @@ -51,6 +61,11 @@ public class GainSchedParameters : ModelParametersBaseClass
/// </summary>
public List<double[]> LinearGains { get; set; } = null;

/// <summary>
/// Threshold for when to use different gains
/// </summary>
public double[] LinearGainThresholds { get; set; } = null;

/// <summary>
/// An array of 95% uncertatinty in the linear gains (u-u0))
/// </summary>
Expand Down
25 changes: 15 additions & 10 deletions TimeSeriesAnalysis.Tests/Tests/PlantSimulatorMISOTests.cs
View file
Original file line number Diff line number Diff line change
Expand Up @@ -127,8 +127,10 @@ public void SetUp()

gainSchedParameters1 = new GainSchedParameters
{
TimeConstant_s = new double[] { 10 },
TimeConstant_s = new double[] { 10,20 },
TimeConstantThresholds = new double[] { 2},
LinearGains = new List<double[]> { new double[] { 5 }, new double[] { 10 } },
LinearGainThresholds = new double[] { 2 },
TimeDelay_s = 0,
Bias = 0
};
Expand Down Expand Up @@ -223,18 +225,21 @@ public void GainSched_Single_RunsAndConverges()
var isOk = plantSim.Simulate(inputData, out TimeSeriesDataSet simData);
Assert.IsTrue(isOk);
SISOTests.CommonAsserts(inputData, simData, plantSim);
double[] simY = simData.GetValues(processModel1.GetID(), SignalType.Output_Y);
double[] simY = simData.GetValues(gainSched1.GetID(), SignalType.Output_Y);

Assert.IsTrue(Math.Abs(simY[0] - (1 * 50 + 0.5 * 50 + 5)) < 0.01);
Assert.IsTrue(Math.Abs(simY.Last() - (1 * 55 + 0.5 * 45 + 5)) < 0.01);
Assert.IsTrue(Math.Abs(simY[N/3-2] - 5) < 0.2,"first step should have a gain of 5");
Assert.IsTrue(Math.Abs(simY.Last() - 30) < 0.2, "third step should have a gain of 10");
// Assert.IsTrue(Math.Abs(simY.Last() - (1 * 55 + 0.5 * 45 + 5)) < 0.01);

/*Plot.FromList(new List<double[]> {
simData.GetValues(processModel1.GetID(),SignalType.Output_Y_sim),
simData.GetValues(processModel1.GetID(),SignalType.External_U,0),
simData.GetValues(processModel1.GetID(),SignalType.External_U,1)
//y[k] = a y[k-1] + b u[k]
Shared.EnablePlots();
Plot.FromList(new List<double[]> {
simData.GetValues(gainSched1.GetID(),SignalType.Output_Y),
inputData.GetValues(gainSched1.GetID(),SignalType.External_U,0),
},
new List<string> { "y1=y_sim1", "y3=u1","y3=u2" },
timeBase_s, "UnitTest_SingleMISO");*/
new List<string> { "y1=y_sim1", "y3=u1" },
timeBase_s, "GainSched_Single");
Shared.DisablePlots();
}


Expand Down
2 changes: 1 addition & 1 deletion TimeSeriesAnalysis.Tests/Tests/PlantSimulatorSISOTests.cs
View file
Original file line number Diff line number Diff line change
Expand Up @@ -236,7 +236,7 @@ public static void CommonAsserts(TimeSeriesDataSet inputData,TimeSeriesDataSet s
Assert.IsTrue(firstTwoValuesDiff < 0.01, "system should start up in steady-state");
Assert.IsTrue(lastTwoValuesDiff < 0.01, "system should end up in stedy-state");
}
Assert.AreEqual(simData.GetLength(), simData.GetTimeStamps().Count(), "number of timestamps shoudl match number of data points in sim");
Assert.AreEqual(simData.GetLength(), simData.GetTimeStamps().Count(), "number of timestamps should match number of data points in sim");
Assert.AreEqual(simData.GetTimeStamps().Last(), inputData.GetTimeStamps().Last(),"datasets should end at same timestamp");

/* foreach (var modelKeyValuePair in plant.GetModels())
Expand Down

0 comments on commit b2e6e00

Please sign in to comment.