Как использовать булевы переключатели в Modelica для предотвращения истощения акций ниже нуля?

Question

Описание проблемы

Я строю библиотеку для поддержки System Dynamics (SD), как моделирование в Modelica. В отличие от свободно доступной библиотеки Cellier et al. Я убежден, что вполне можно использовать акаузальные соединители: передача стоимости запаса в качестве "потенциала" через соединители позволяет создавать компактные компоненты (например, потоки = процессы).

В SD мы могли бы отличить материальные ("массовые") запасы от информационных, которые могут стать отрицательными. Чтобы поддержать это, я использую следующие определения для массового порта (здесь дается определение для StockPort - его аналог, FlowPort, просто будет иметь логические входные переменные вместо выходных переменных и будет дано позже):

 connector StockPort "Used to represent stock and flow connections"
    Real stock "Current value of material in the stock";
    flow Real rate "Flow that affects the stock";
    // Boolean switches
    output Boolean stopInflow "True indicates that nothing can flow into the stock";
    output Boolean stopOutflow "True indicates that nothing can flow out of the stock";
 end StockPort;

Булевы переключатели указывают для каждого порта запаса, разрешено ли заполнение или слив.

Для " материального запаса " stopOutflow Переключатель должен предотвращать слив запаса ниже нуля. К сожалению, в следующем примере это не сработает: запас опустится чуть ниже нуля.

Минимальный пример использования разъемов

Следующие TestModel использует эти строительные блоки:

• function constrainedRate( indicated rate, stopInflow, stopOutflow) используется для возврата скорости, которая будет соответствовать заданным ограничениям (т. е. булевых переключателей)

• connector StockPort как описано выше

• connector FlowPort аналог к StockPort
• model MaterialStock складской компонент с одним StockPort который не должен быть дренируемым ниже нуля
• model LinearDecline элемент потока с одним FlowPort (т. е. Sink), который моделирует слив подключенного запаса с постоянной скоростью (здесь установлено значение 1)

Основная модель просто инициирует stock с initialValue = 5 это связано с process линейного спада с declineRate = 1,

model TestModel "Stop draining a stock below zero"

  function constrainedRate "Set rate for a port according to signals from stock" 
    input Real indicatedRate "Proposed rate for port of flow element";
    input Boolean stopInflow "Signal from connected stock";
    input Boolean stopOutflow "Signal from connected stock";
    output Real actualRate "The rate to use";
  protected
    // check whether indicated rate is negative (e.g. an inflow to the connected stock)
    Boolean indRateIsInflow = indicatedRate < 0;
  algorithm
    // set rate to zero if stopSignal matches character of flow
    actualRate := if indRateIsInflow and stopInflow 
          then 0 
        elseif not indRateIsInflow and stopOutflow 
          then 0 
        else indicatedRate;
  end constrainedRate;

  connector FlowPort "Used to represent stock and flow connections"
    Real stock "The current stock level (e.g. Potential) of a connected stock or flow data for special stocks";
    flow Real rate "Flows that affect the material stock";
    input Boolean stopInflow "True indicates that nothing can flow into the stock";
    input Boolean stopOutflow "True indicates that nothing can flow out of the stock";
  end FlowPort;

  connector StockPort "Used to represent stock and flow connections"
    Real stock "Current value of stock";
    flow Real rate "Flow that affects the stock";
    output Boolean stopInflow "True indicates that nothing can flow into the stock";
    output Boolean stopOutflow "True indicates that nothing can flow out of the stock";
  end StockPort;

  model MaterialStock "Stock that cannot be drained below zero"
    StockPort outflow;
    parameter Real initialValue;
  protected
    Real x(start = initialValue);
  equation
    // rate of change for the stock
    der(x) = outflow.rate;
    // ports shall have level information for stock
    outflow.stock = x;
    // inflow to stock is unrestricted
    outflow.stopInflow = false;
    // provide Boolean signal in case of negative stock
    outflow.stopOutflow = x <= 0;
  end MaterialStock;

  model LinearDecline "Decline of stock at a constant rate"
    FlowPort massPort;
    parameter Real declineRate(min = 0) "Rate of decline (positive rate diminishes stock)";
  protected
    // a positive rate should drain the stock (which here matches Modelica's rule for flows)
    Real rate(min = 0);
  equation
    rate = declineRate;
    // observe stock signals and constraints
    assert(rate >= 0, "Rate must be positive and will be set to zero", level = AssertionLevel.warning);
  // set the rate according to constraints given by stock
    massPort.rate = constrainedRate( max(rate, 0), massPort.stopInflow, massPort.stopOutflow );
  end LinearDecline;

  // main model
  MaterialStock stock( initialValue = 5 );
  LinearDecline process( declineRate = 1 );
equation
  connect( stock.outflow, process.massPort );
end TestModel;

Имитация модели с использованием DASSL из StartTime = 0 в StopTime = 10 показывает ожидаемое поведение для переменной stock.outflow.stock:

К сожалению, значение немного ниже нуля при t = 5.0 и позже.

Каким-то образом событие (стоимость акций <= 0) обнаруживается слишком поздно. Что я могу сделать?

(До сих пор Imo неликвидное лекарство было использовать when событие (state-event) для reinit стоимость акций до нуля. Мои эксперименты с использованием noEvent обертки на if операторы и логические условия также не были успешными.)

Answer 1

Не могли бы вы проверить это решение? Использование Modelica.Constants.eps у меня работает. Я также использовал Dassl, и он работал для разных размеров. Я изменил следующую строку с:

// provide Boolean signal in case of negative stock
outflow.stopOutflow = x <= 0;

в

// provide Boolean signal in case of negative stock
outflow.stopOutflow = x <= Modelica.Constants.eps;

Затем выходной массив (рассматриваемый в python для этого поста):

Output using a zero instead of eps:
[ 5.00000000e+00  4.00000000e+00  3.00000000e+00  2.00000000e+00
  1.00000000e+00  0.00000000e+00 -7.37276906e-12 -7.37276906e-12
 -7.37276906e-12 -7.37276906e-12 -7.37276906e-12 -7.37276906e-12
 -7.37276906e-12 -7.37276906e-12]
Output using eps:
[5. 4. 3. 2. 1. 0. 0. 0. 0. 0. 0. 0. 0.]