=Paper=
{{Paper
|id=Vol-2289/paper18
|storemode=property
|title=Diagnosing Hybrid Cyber-Physical Systems using State-Space Models and Satisfiability Modulo Theory
|pdfUrl=https://ceur-ws.org/Vol-2289/paper18.pdf
|volume=Vol-2289
|authors=Alexander Diedrich,Oliver Niggemann
|dblpUrl=https://dblp.org/rec/conf/safeprocess/DiedrichN18
}}
==Diagnosing Hybrid Cyber-Physical Systems using State-Space Models and Satisfiability Modulo Theory==
Diagnosing Hybrid Cyber-Physical Systems using State-Space Models and
Satisfiability Modulo Theory
Alexander Diedrich1 and Oliver Niggemann2
1
Fraunhofer IOSB-INA, Lemgo, Germany
2
Institute Industrial IT, Lemgo, Germany
e-mail: alexander.diedrich@iosb-ina.fraunhofer.de
oliver.niggemann@hs-owl.de
Abstract alarms can quickly overwhelm an operator and thus keep
him from finding the faulty component that caused the fault
Diagnosing faults in large cyber-physical produc- [2].
tion systems is hard and often done manually. Identifying and isolating faults in large industrial plants can
In this paper we present an approach to lever- take precious time which can quickly lead to a significant
age methods from the fault detection and isolation deterioration of the produced material or even physical de-
community as well as model-based diagnosis to struction of involved components. For operators of these
diagnose faults. Given a model of the production plants this can lead to high costs. These high costs occur
system we capture its dynamic behaviour with a even in smaller scale enterprises, due to low quality output,
state-space model. Then we translate the state- costs to locate and repair the broken component, and costs
space model into satisfiability theory modulo lin- to ramp up production again after the fault.
ear arithmetic over reals. This translation converts In the presented approach we attempt to perform fault de-
numerical information in symbolic logic. These tection and isolation (FDI) through model-based diagnosis
symbols can be used to diagnose faults with Re- over satisfiability modulo theory (SMT). For this, a model of
iter’s diagnosis algorithm. the physical process is created manually. Industrial cyber-
We use a four-tank model as a demonstration use physical production systems are dynamic because their pa-
case. Under the assumption that the use-case is rameters change over time and contain hybrid signals which
fully-observable (i.e. all components except the can be either binary or continuous. Here we use state-space
water tanks can be observed) our methodology de- models to capture the dynamic behaviour of hybrid and
tects all injected faults. cyber-physical production systems. We limit our use case
to a multiple tanks model. The behaviour of the tanks is
1 Introduction modelled using differential equations, while the connection
between components is modelled using predicate logic. A
For operators of large industrial cyber-physical production
set of piecewise functions translates the state space model
systems (CPPS) it is often a hard task to precisely detect,
into satisfiability theory modulo linear arithmetic over re-
identify, and isolate technical faults [1]. This is especially
als (LRA). Through the translations it is possible to lever-
the case in large production plants in the process industry or
age standard model-based diagnosis algorithms, once devel-
in pharmacological processes which often extend over sig-
oped to diagnose binary circuits, to diagnose hybrid cyber-
nificant physical distances, consist of highly interdependent
physical production systems. The outcome of this trans-
components, and involve many parallel paths in the form of
lation is a set of tuples which states for each component
pipe and valve networks. The correct behaviour of a biolog-
whether or not it is currently faulty < comp, status >.
ical reactor, for example, depends on the exact amount of
This tuple is converted into predicate logic and combined
different ingredients, their pressure, their temperature, their
with the connection model of the plant. We employ Reiter’s
viscosity, the ambient temperature and pressure, or other im-
diagnosis lattice [3] to find the minimum cardinality diagno-
portant process values. Further, subsequent processes heav-
sis and thus isolate the fault that caused the production plant
ily depend on the correct amount, quality, and time of dis-
to fail.
charge from the biological reactor. Since these systems are
In this work we demonstrate how our described approach
interdependent it is hard for human operators to physically
can be used with a multiple tanks model, how we translate
locate the root-cause of a fault. For example, a valve might
the numerical state-space model into predicate logic, and
break and block the flow of some liquid into the reactor. Due
how we can perform diagnosis using Reiter’s well-known
to this blocked flow, the pressure and temperature within the
algorithm [3]. We show that in case of a fully-observable
reactor might change and degrade the material by changing
system our approach is able to find all faults as part of
its viscosity. The degraded material might then go on into
the minimum cardinality diagnosis. By fully-observable we
a stamping and forming process. Consequently, due to the
mean that we can observe the behaviour of all components
changed viscosity, the presses will register changed pres-
except the water level within the tank. The water level is
sures within their control systems. For a human operator,
inferred through calculation in the state-space model. We
all these components will at some point sound an alarm in-
also give ideas how the approach can be extended to semi-
dicating that parameters are outside their normal operating
observable and non-observable systems.
conditions. In modern industrial plants, the amount of these
This paper makes the following contributions: Assumption 2 (Measurement errors). Measurements from
1. We show how to capture timing behaviour for model- the flow sensors and over- and underflow switches are al-
based diagnosis on the basis of hybrid and dynamic ways perfect without measurement error. If necessary, it is
cyber-physical systems. stated when this assumption is relaxed.
2. We demonstrate how diagnosis methods (i.e. Reiter’s Justification for assumption 1 is that it suffices to simu-
algorithm) from the model-based community can be late faults within the valves and tanks. Modelling the pipes
successfully combined with approaches from the fault with physical properties would dramatically increase the
detection and isolation (FDI) community. model size and thus reduce the clarity. Assumption 2 is
taken to simplify the used equations. When necessary this
3. We show how satisfiability modulo theory can help to assumption is relaxed.
perform diagnosis in hybrid and dynamic systems by This demonstration use case can be imagined as a prepro-
keeping the amount of computations low. cessing stage in a larger industrial plant within the process
industry. A reliable external water supply is provided by
2 Demonstration Use Case the facilities of the industrial park. The water flows from
the supply line into a buffer tank t0 . From there it can go
For this work we will use the four tank system depicted in
into one or both or the two intermediate tanks or bypass
Figure 1 as a running example. The system consists of four
both tanks to go directly into tank t3 . The two intermediate
tanks can be thought of as a mixing stage (not modelled)
where ingredients are added to the water until it reaches the
holding tank t3 . From the holding tank the water flows to
subsequent process steps.
p0
The water level in tanks can be described by well-known
t0 differential equations. Laubwald [4] provided a comprehen-
sive overview about modelling multiple water-tank systems.
A single tank can be described with the differential equation
p1 p2 dh
Qi (t) − Qo (t) = A (1)
dt
which describes the time derivative of the height h given
p3 t1 t2 the tank area A. Qi is the inflow to the tank and Qo is the
outflow. However, in the real world tanks are subjected to
gravity and properties of their materials. Therefore the out-
p4 p5 flow of a tank is calculated by
√
Qo = Cd a 2gh (2)
t3 where Cd is the discharge coefficient taking into account all
fluid characteristics, losses, and irregularities and g is the
p6 gravitational constant. a is the cross sectional area of the
orifice within the tank. All tanks have a perfectly circular
bottom and a cylindrical shape. Combining equations 1 and
Figure 1: The demonstration use case showing a four-tanks 2 one can create
model √ dh
Qi (t) = Cd a 2gh + A (3)
water tanks t, seven electric valves p with integrated flow dt
sensors, an unlimited water source and an unlimited water and reformulating this to bring dh
dt on the left-hand side
sink (not shown). Valve p0 controls water from the unlim-
dh 1 √
ited water source, for example the public water mains, into = (Qi − Cd a 2gh) (4)
tank t0 . From there three pipes with an equal diameter lead dt A
to valves p1 , p2 , and p3 . Valve p1 leads into tank t1 and To calculate a new tank height h at time t, given the previous
valve p2 leads into tank t2 . Valve p3 bypasses both tanks height h0 one can use
and is directly connected to tank t3 . Over valves p4 and p5
the two tanks can be drained into tank t3 . Finally, valve p6 h(t) = h0 + ∆h (5)
drains tank t3 into the unlimited water sink, for example a Through substitution into equation 4 this results in
river or a processing facility.
Each tank has two binary sensors which indicate overflow 1 √
h(t) = h(t − 1) + (Qi (t) − Cd a 2gh(t − 1)) (6)
and underflow, respectively. There are no provisions to di- A
rectly measure the water level. Each valve has a switch which describes that, given the discharge coefficient, the
which indicates whether or not the valve is enabled. In ad- gravitation, and the diameter of the orifice it is possible
dition, each valve has an associated flow sensor. For the to calculate a new height from a given input with only the
present system the following assumptions are made: knowledge of the previous water level.
Assumption 1 (Pipes). Pipes are invisible to the system, These differential equations can be used to create a simu-
have ideal physical properties and cannot break lation of the four-tank system depicted in 1. The following
Parameter Element Value Unit and Mealy automata.
Area Daigle et al. [7] have adapted a discrete event approach to
A t0 4.0 m2 diagnose continuous systems. They claim that each fault
t1 2.0 m2 that occurs in a continuous system has a unique fault sig-
t2 2.0 m2 nature. A fault signature denotes a qualitative effect that a
t3 6.0 m2 fault occurs in an observation. They also claim that there ex-
Height ists a measurement ordering that describes which sequence
H t0 20.0 m measurements deviate until a fault occurs. To capture fault
t1 10.0 m ordering they manually construct a temporal causal graph.
t2 10.0 m Under the assumption that all fault signatures and measure-
t3 20.0 m ment orderings are known, they employ a diagnoser that
Discharge Coefficient traces the states through the temporal causal graph based on
Cd t0 1.0 None measurements. The output of the diagnoser is a fault trace.
t1 1.0 None A second diagnosis algorithm takes this fault trace and de-
t2 1.0 None termines which components must be faulty to explain this
t3 1.0 None trace. This second diagnosis step is similar to the diagnosis
Orifice lattice introduced by Reiter.
a p0 0.3 m2 Grastien et al. [8] have developed an approach to extend
p1 0.1 m2 Reiter’s diagnosis algorithms which was described for bi-
p2 0.1 m2 nary circuits to include discrete event systems and hybrid
p3 0.1 m2 systems. Their approach is similar to Daigle et al, Struss,
p4 0.1 m2 and Provan in that they transform the continuous parts of a
p5 0.1 m2 model into qualitative states. Following this, their preferred-
p6 0.3 m2 first algorithm goes through Reiter’s diagnosis lattice and
Gravitational Constant g 9.81 m/s
2 computes valid hypothesis with the goal of finding a mini-
mum cardinality diagnosis. An improvement over previous
Table 1: Parameters of the four-tank system for the demon- work is that they implement their hypotheses tests with a
stration use case SAT solver.
Roychoudhury et al. [9] [10] have shown how to use hybrid
bond graphs (HBG) to diagnose hybrid systems. HBGs ab-
parameters were used to run the system simulation: The pa-
stractly model the system by describing causal, continuous
rameters in table 1 have been chosen to represent a typical
relationships between components. In Daigle et al. [7] they
industrial use case. At the same time, noise parameters such
have shown how to employ the developed HBGs to diag-
as the discharge coefficient are kept neutral to further the ar-
nose a spacecraft power distribution system. Prakash et al.
gument. The parameters area and height for each tank have [11] have used an extended framework with HBGs to make
been chosen such that they are big enough to hold and store
improvements in diagnosing two-tank systems.
some water while the experiments are running. This way
Grastien [12] used SMT for the diagnosis of hybrid systems.
occurring errors will also have a longer time to propagate.
He discretizes values in a hybrid system into a set of dis-
Further, the orifices for each tank are quite small. This pro-
tinct states. Each observation < τ, A >is understood as a
longs the time it takes to drain the tanks in case the water
behaviour A at time τ , where A is a partial assignment of
supply stop through the occurrence of a fault.
the variables in a state. Measurement errors are included by
including constraints which state that the observed voltage
3 Related Work must be between two tolerance thresholds. Each variable is
Struss [5] published a paper on the fundamentals of MBD of augmented with an indicator stating at which time-step the
dynamic systems. In this he described how hybrid systems variable expression is valid.
can be modelled without resorting to a complete simulation Fraenzle et al. [13] have augmented SMT with stochas-
of the system under investigation. He proposed to capture tic in order to analyse stochastic hybrid systems. By using
the temporal and dynamic behaviour of a hybrid system in a bounded-model checking together with probabilistic hybrid
set of modes which model the system. Each mode has dis- automata, piecewise deterministic Markov processes, and
tinct state and temporal constraints in addition to so called stochastic differential equations they are able to create an
Continuity, Integration, and Derivatives (CID) constraints analysis system without the need to formulate intermediate
that affect all modes. For one mode, all variables have a do- finite-state abstractions as the methods mentioned above do.
main which captures the permissible states (i.e. values) for In another work Khorasgani and Biswas [14] describe a
this variable. Diagnosis is performed by checking whether hybrid system model through hybrid minimal structurally
the set of constraints together with the observations from overdetermined sets (HMSOs). These are sets of differential
sensors is consistent. He demonstrates his approach on a equations and (in-) equations which model the behaviour of
car’s anti-braking system and claims to find all usually oc- a hybrid system. Their FDI algorithm works as follows: The
curring faults. algorithm detects the current system mode and generates an
When dealing with hybrid systems there always exists the appropriate model. From this it generates a minimal set of
problem of discretization. Provan used a composite au- HMSOs for this mode. The residuals are computed for each
tomaton for this. Struss divides his system into modes by HMSO and can then be combined with fault signature to
discretizing the underlying sensor values. Lin [6] already perform diagnosis.
showed in 1994 that online and offline diagnosis for discrete In contrast to Struss, Provan, and Lin we do not use au-
event systems (DES) can be realised by using simple Moore
tomatons and mode estimation to partition the system into 8. To perform model-based diagnosis according to the prin-
different states. Instead, we only sample the system at some ciples proposed by Reiter [3] it is necessary to separate the
suitable interval and use the obtained information directly diagnostics part from the state propagation. Therefore, we
to model the states in the state-space representation. Un- will not calculate classical residuals as in equation 8 and can
like to space-craft in the case of Daigle in industrial systems further remove the gain factor and residuals Kr(t).
fault signatures and measurement orderings are unknown, For this work the state-space model needs to be described
which requires us to pursue a more uninformed approach. more abstractly. First the top-level information flow is de-
Our approach can be more seen as an alternative to hybrid scribed. This shows how the state is propagated within the
bond graphs used by Roychoudhury et al., while they are system. Here, the model is general enough to be extended
at the same time an extension to the work of Grastien and and adapted to many use cases. After that the state-space
Khorasgani and Biswas. In comparison to Grastien we do model is described on the demonstration use case introduced
not singly use satisfiability modulo theory, but instead cap- in section 2. In a third step the diagnosis part is described
ture system behaviour in a state-space representation. We which fusions the calculation of binary residuals with an ex-
expect this to reduce the required computational effort. We pression in predicate and SMT logic. The state is propagated
also make use of (in-) equations and differential equations through
as were used by Khorasgani and Biswas, but augment these
x(t + 1) = f (x(t), u(t))
with the diagnostic reasoning of traditional model-based di- (9)
agnosis. Compared to Fraenzle, we do not make use of y(t) = g(x(t), u(t), τ )
stochastic SMT at this point to keep the system more ex- where x(t + 1) is a vector of the state in the next time step,
plainable for users. x(t) is the current state vector, u(t) is the observable input
vector, y(t) is the observable output vector, and τ is a vector
4 Modelling a Hybrid System of threshold values.
According to the demonstration use case the water level can-
This section first shows the requirements for developing and
not be measured directly. Therefore, each tank’s water level
evaluating a FDI method for hybrid, cyber-physical systems.
needs to be calculated through its inflow and outflow. The
Then it shows the concept to realise these requirements. The
inflow and outflow can be measured through the associated
developed approach makes use of MBD by modelling the
valves in each in- and outflow pipe. Since each tank has
hybrid system with a state-space representation. This model
sensors to indicate under- and overflow, these are used for
is augmented with an observer, which determines boolean
the target (output). For the state, input, and output vectors
residuals. These residuals indicate whether or not a compo-
we thus have:
nent is faulty. The information about which components in-
dicate faults are merged with sensor observations and make overf low0
up the diagnostic part of the approach. Diagnosis is done ..
h0 f low0 .
through Reiter’s diagnosis lattice. f low1
h1 overf low
The formal form of a state-space representation is: x= u= .. . y =
3
h2 underf low0
h3 ..
x(t + 1) = (A + ∆A)x(t) + (B + ∆B)u(t)+ f low6 .
Bd dk (t) + Bafa (t) + Bc fc (t) (7) underf low3
y(t) = (C + ∆C)x(t) + Ds fs (k) + Dω ω(k)
The function f (x, u) models the current state and its current
x(t) is the systems’s state, u(t) the control input, y(t) input and from this computes the next state. Therefore we
the observed output, fa (t) the unexpected actuator fault, can write:
fc (t) a component fault, fs (t) a sensor fault dk (t)
a process disturbance, and ω(t) measurement noises. f (x(t), u(t)) = A∆(x, u, t) + Bu(t) (10)
A, B, C, Bd , Ba , Bc , D, s, Dω are known parameter matri- with the connection matrices being
ces and ∆A, ∆B, ∆C modelling parameter errors.
For observer-based methods and to calculate residuals these 1 0 0 0
0 1 0 0
equations can reformulated to: A=
0 0 1 0
x̂(t + 1) = Ax̂(t) + Bu(t) + Kr(t) 0 0 0 1
r(t) = y(t) − ŷ(t) (8) and
ŷ(t) = C x̂(t)
1 −1 −1 −1 0 0 0
0 1 0 0 −1 0 0
B=
0
Here, x̂(t) and ŷ(t) are estimates of the state and output val-
0 0 1 0 0 −1
ues. r(t) is the calculated residual signal, and K is a gain
0 0 0 1 1 1 −1
factor.
According to assumption 2 we can safely neglect the factors A shows that each current state only influences the exact
the ∆ and ω terms in equation 7. Further, in this approach next state. For the demonstration use case this means each
we will only model component faults. Therefore, we can differential equation which models a tank will only affect
remove fa (t), fs (t), and through assuming that there are no the state of this single tank. Matrix B shows the connec-
disturbances within the process we can also remove dk (t). tions between the system’s components, which in this case
This leaves only the system’s observable input, observable are the pipes between the tanks. The first row describes how
output, state, and component fault in equation 7. Through the system’s primary input is connected as an input (indi-
these simplifications equation 7 becomes closer to equation cated by the number 1 in the first column) to tank 1. The
three values of −1 in the first row show there are three pipes components. In the presented demonstration use case two
that are used as the output of tank 1. fault models for tanks and valves exist, respectively. Tanks
To model the water level in each tank it is possible to use dif- fail, when the water level within the tank reaches either
ferential equations. Each differential equation has the form the upper limit (overflow) or the lower limit (underflow).
introduced in section 2: Pumps fail when the measured flow deviates more than a
1 √ certain amount from the expected flow.
h(t) = h(t − 1) + (Qi (t) − Cd a 2gh(t − 1)) (11) Classical MBD uses observations(OBS), a system
A
description(SD), and a component description(COMPS)
Using the parameters from the state-space system this is for describing a system. After having described the
written as actual behaviour of the hybrid system with state-space
1 √ equations, it is now important to translate this into diag-
x(t + 1) = (u(t) − Cd a 2gx(t)) (12)
A nostic information. OBS are given by the input vector
A vector ∆(x, u, t) can be created with the right-hand side u(t). The component behaviour COMPS is described
of these equations: by the differential equations in the case of tanks and by
√ assuming no further properties for the valves, resulting
x0 (t + 1) = A10 (u0 (t) − Cd,0 a0 2gx0 (t)) in input(valvei ) = output(valvei ). SD is given in two
x (t + 1) = 1 (u (t) − C a √2gx (t)) parts. For normal operation this is the incidence matrix B, a
1 1 d,1 1
√ 1
∆(x, u, t) = A1
predicate logic description of the inputs and outputs of the
x2 (t + 1) = A12 (u2 (t) − Cd,2 a2 2gx2 (t))
√ system, and a fault model. For the given demonstration use
x3 (t + 1) = A13 (u3 (t) − Cd,3 a3 2gx3 (t)) case it suffices to specify a weak-fault model (WFM). A
With this model it is possible to propagate the state of WFM to model the fault modes of the tanks can be specified
the system as it evolves through time. Differential equa- as
tions calculate the water level in the tank for the next state, σT,i : HT,i → ¬o(i) ∧ ¬l(i) (16)
given the current water level and the inflow obtained by For valves the statement is specified as
reading the valve flow sensors. However, given this infor-
mation a control system cannot yet determine the full be- σP,i : HP,i → (f lowil ≤ f lowi ) ∧ (f lowiu ≥ f lowi ) (17)
haviour of the system. For this, the output vector y(t) = In this case the health variables H do not describe a proba-
g(x(t), u(t), τ ) needs to be calculated. bility for the component being faulty, but are instead binary.
The terms σt,i and σp,i can be written as a vector
o(h0 , τ0o )
.. T
. C = [σT,0 ... σT,3 σP,0 ... σP,6 ]
o(h , τ o )
3 3
g(x(t), u(t), τ ) = C (13) If C is semantically interpreted through an SMT solver C ′ =
l(h0 , τ0l )
.. I(C), we obtain the diagnosis vector
.
C′ = [c0
T
l(h3 , τ3l ) c1 ... c10 ]
and with ci ∈ {⊤, ⊥}. This vector shows for each component
1 0 0 0 0 0 0 0 whether it is faulty or not, given the current observations
0 1 0 0 0 0 0 0 from the sensors.
0 0 1 0 0 0 0 0
The numerical information for the statements σT,i and σP,i
0 0 0 1 0 0 0 0
C= is obtained from the state space model. Within the state-
0 0 0 0 1 0 0 0 ments the state vector x(t) represents the water level hi and
0 0 0 0 0 1 0 0 the input vector u(t) represents the flow values f lowi . By
0 0 0 0 0 0 1 0 interpreting the statements it is possible to translate the sub-
0 0 0 0 0 0 0 1 symbolic, numerical data within the state-space model into
τ is a vector of threshold values which indicate at what symbolic information used for diagnosis through the vector
height the tank is overfull or underfull. For notation we use C′
τio to denote the threshold for the upper limit of tank i and Equation 9 shows the propagations of the state vector
τil to denote the lower limit of tank i. The diagonal matrix through time. For each new time step the statements 16
C maps the results of the functions o(h, τ ) and l(h, τ ) into and 17 have to be reformulated. To capture this time-related
the output vector y. The function o(h, τ ) indicates when the behaviour we adopt the notation of Grastien [12] and state
water level within the tank has approached the upper limit. that varname@t stands for the variable varname at time
This is calculated by t, where t ∈ N0 . From this, we can state the value for each
{ variable at each observed time step. When the observations
0 if hi ≤ τio are only carried out while the system is in normal operation
o(hi , τio ) = (14) it is possible to create a logical representation of all obser-
1 else
vations so far:
Likewise, the lower limit of the water level can be calculated ∧ ∧
{ σT,i @t ∪ σP,i @t
0 if hi ≥ τil
l(hi , τil ) = (15) t t
1 else
which describes the logical conjunction of all σT,i and σP,i
To diagnose faults within the described state-space system over all time steps. Adding the statements in each time
it is necessary to obtain health information about single step to the knowledge base as described by Grastien will
increase the required space linearly and still take exponen- that we can measure the water flow through each valve
tial time to check the consistency. Especially in large in- at each point in time. This is also a realistic assumption
dustrial plants where observations run for months with in- for many smaller scale application and most demonstration
dividual observations being performed at second intervals, plants which are built with observability in mind. Older and
a linearly growing knowledge base is infeasible. For exam- more complex industrial plants, however, contain more of-
ple, when observing 200 signals with a sampling rate of one ten component whose parameters can not be observed.
second a knowledge base would grow by 17,280,000 data
points per day. Therefore, in this work we will focus only Semi-observable system
on the observations in the current time step. This keeps the In the case that the system is semi-observable, not all com-
knowledge base size constant and adds no additional com- ponents’ behaviour can be observed with sensors. This is the
putational complexity. More observations can bring a higher case in most real-world industrial plants where its either too
precision in locating a fault. In this case, the number of ob- expensive to add sensors for every machine parameter or it
servations can be increased by some constant factor, taking is infeasible due to physical constraints or historical reasons.
for example always the last three observations into account. A diagnostic system which cannot observe every parameter
A hybrid system can be represented through a directed- has to work with partial information. If necessary the miss-
acyclic graph (DAG) showing the connections and causal ing values need to be estimated. In the case boolean circuits
relationships between components. Depending on the loca- this can be done straightforward. For every component in a
tion of the component within the graph a fault in one com- boolean circuit the behaviour model is known. Further, the
ponent may cause several other components to fail as well. system description SD is known. If only parts of the com-
In the demonstration use case, for example, if valve p0 fails, ponents can be observed a diagnostic reasoning system can
all the other valves and the tanks will also exhibit anoma- infer the missing values.
lous behaviour. The goal in diagnosis is therefore to find In hybrid, physical systems inferring values is more dif-
the smallest amount of components which would explain a ficult. Some real-world components may behave non-
fault. This search for minimum cardinality diagnoses can linearly, stochastically, or very unpredictable. Further, sig-
be done with Reiter’s diagnosis lattice. First, a power set nal propagation may not be instantaneous as in boolean cir-
P = P(COMPS) is constructed. This contains all sets of cuits, but a change in one parameter may only be notice-
sets of components. From this the diagnosis lattice can be able some time later. This is the case, for example, in bio-
created. On the bottom is the empty set which denotes no reactors. If the temperature in a reactor changes, the sub-
faulty components. In the row above are all sets that con- stance may only exhibit a change in an observable property
tain exactly one component. In the row above that are sets some time later.
that contain exactly two components and so forth. Each In the demonstration use case the tanks are modelled
observation of the sensors within the system leads to a re- through differential equations and the valves have the
computation of the set of possible faulty components C ′ . By throughput that is maximally allowed by the outflow of a
computing the hitting sets of all these observations it is pos- tank. Thus, even if not all sensors can be observed it is still
sible to close in on the faulty components. In the diagnosis possible to infer missing values.
lattice this is done by searching the lattice bottom-up and Non-observable system
refuting all branches which include a component that can be
In non-observable systems only the primary inputs and out-
proven to be healthy through observations. Once the lattice
puts are observable. A diagnostic system needs to measure
has been searched the solutions with the minimum number
the primary input signals, the primary output signals and
of components are the minimum cardinality diagnoses ω ′ .
combine these with the model SD and COMPS of the hy-
Once the diagnosis framework has been set up three possible
brid system. To perform diagnosis, every intermediate value
usages can be identified: fully-observable, semi-observable,
must be assumed by propagating the primary input values
and non-observable. In the first type of usage the diagnoser
through the system. This approach is the most computing
can observe every property of the physical system, except
intensive, since assumable values need to be computed in se-
the water level described by the state vector x(t). In the sec-
quence. Diagnosis is performed by comparing the expected
ond type of usage, only a subset of sensors are accessible to
primary output values with the observed primary output val-
the diagnoser. Thus, some values need to be approximated.
ues and then going through the circuit back-to-front to find
In the third type of usage only the primary inputs and pri-
diagnosis candidates.
mary outputs can be observed, while all other values need to
be inferred. The following three sub-sections describe these
types of usage in detail. 6 Experiments
To show that the developed diagnostic methodology works
5 Observability as intended 16 experiments with the simulation of the
In this paper we will focus on diagnosis of fully-observable demonstration use case were carried out. These are divided
systems. This simplifying use case makes it easier to de- into two sets. In the first set the primary input to the demon-
scribe the methods, while the extension to semi-observable stration use case was a constant stream of water. We expect
systems, and non-observable systems is reserved for future in this case that during normal operation the water level in
work. However, this section gives some ideas on how to ex- the tanks will reach a constant height and remain there until
tend the developed methodology to include semi- and non- a fault occurs. In the second set of experiments the primary
observable systems. input was changed to a sinusoidal water stream. For this we
used the function
Fully-observable system {
In full-observable systems we assume that each component O + α sin(2t) if t ≤ T π
in(t) = (18)
can be observed. For the demonstration use case this means 0 else
Equation 18 shows the form of the sinusoidal wave with pe- Constant Sinusoidal
riod T . We use an offset C to ensure a constant basic input Index # Faults Detected # Faults Detected
stream into tank 1. The gain factor A adjusts the period such 0 p0 x* p0 x*
that we achieve variability within the tank water levels, but 1 p3 x* p3 x*
without triggering and under- or overflow. Further, we use 2 p5 x* p5 x*
a piecewise function to cut off the negative half-wave of the 3 p6 x* p6 x*
sinusoidal input wave for convenience and ease of interpre- 4 p1 , p3 x* p1 , p3 x*
tation. 5 p4 , p5 x* p4 , p5 x*
For both sets of experiments the normal operating condition,
five cases of single-faults, and two cases of double-faults Table 2: Recognized faults for experiments with constant
were simulated. For each experiment 300 time steps were input stream or sinosoidal input stream at time-step 100
carried out, with the respective faults being injected at time
step 100 and being removed at time step 200. Constant Sinusoidal
We split the developed method into two parts. Part one is Index # Faults Detected # Faults Detected
the quantitative simulation of the demonstration use case 0 p0 x (11) p0 x (11)
described in section 2. Part two is the diagnosis algorithm 1 p3 x (3) p3 x (3)
consisting of the state-space model, the SMT logic, and the 2 p5 x (3) p5 x (3)
diagnosis lattice. Both parts have been written in Python 3 p6 x (3) p6 x (3)
3.4.5. The quantitative simulation provides the user with 4 p1 , p3 x (5) p1 , p3 x (7)
functionalities to inject faults and generate normal process 5 p4 , p5 x (6) p4 , p5 x (6)
data. The location and number faults can be specified as
well as the type of input (for example, if the water inflow Table 3: Recognized faults for experiments with constant
is constant or sinusoidal). The output of the simulation is a input stream or sinosoidal input stream at time-step 199
.csv file which contains all process data as well as the diag-
nostic information. This method was chosen to be close to
real industrial use cases. if valve 5 stops working its flow would almost immediately
SD is modelled through predicate logic. In the fully- go to 0. The sampling frequency is high enough to detect
observable case for the demonstration use case it suffices to this decrease in the flow rate early enough that the water
explicitly model the connections between components, in- level in the tanks is not yet significantly affected. However,
puts, and outputs. Therefore, only three functions are used in large industrial plants sampling rates are often far lower.
for the predicate logic: The function component(c) with A faulty component might then only be recognised once its
arity 1, and the function input(i, c) and output(o, c) with effects have propagated into other observations from other
arity 2. These model the names of components in the sys- components. Further, in the semi- and non-observable cases
tem and the number and names of their inputs and outputs. not every status of every component is known. In this case,
In addition the relation connects(ci , cj ) specifies which in- too, the set of possible faulty components will grow in size.
put is connected to which output. For the present use case As the criterion in table 2 we evaluated the output of the di-
the constants source and sink are used to denote primary agnosis algorithm in the time step 101 which was directly
inputs and primary outputs, respectively. With this logic it after the fault had been injected. It is evident that due to
is possible to describe the connections between components the SMT logic statement in equation 17 every unexpected
in the form: change in the throughput of a valve would be immediately
component(t0)∧ detected.
However, table 3 shows the results when the output of the
... diagnosis algorithm was evaluated directly before removing
component(p6)∧ a fault at time-step 199. The number in brackets denotes the
input(t0.i0, t0) size of the minimum-cardinality set, while x still denotes
whether or not the fault was within the result set. As was
...
the case in table 2, all faults were detected, though the re-
output(p6.o0, p6)∧ sults set also contained components which were not faulty.
connects(source, p0.i0) Figure 2 shows the development of the water level in tank
... 3 for the experiment with a constant water input and the
connects(t3.o0, p6.i0) fault being injected at valve p5 . In this case, the flow from
tank 2 into tank 3 is blocked. In the figure the grey line
connects(p6.o0, sink) represents normal working behaviour and the black line ab-
normal working behaviour. The shaded area indicates the
7 Results time during which p5 is simulated to be faulty. From time-
Table 2 shows the experiments for constant and sinusoidal step 0 until 100 it is evident that the water level in tank 3 is
input streams, the injected fault and whether or not the fault the same in both conditions. The tank was initialised with
was detected. An x in the column detected denotes that the a water level of 7m, which first drains as not enough water
injected fault was among the result set of the diagnosis algo- is flowing into the tank. Then, once tanks 0, 1, and 2 have
rithm. This means the algorithm is complete. An x∗ denotes reached their normal operating conditions tank 3 stabilises.
that only the injected fault was detected, which corresponds Once the fault is injected the water level in tank 3 becomes
to soundness of the algorithm. It must be noted here, how- unsteady. This results from the implementation of the model
ever, that finding only the injected faults depends heavily on which only models tank levels until they reach their upper
the granularity of the underlying data source. For example, limit. Thus, with p5 being disabled tank 2 begins to over-
flow. This in turn redistributes the water pressures in all the [2] Marta Fullen, Peter Schüller, and Oliver Niggemann.
other tanks. Once the fault is removed, however, all flows Defining and validating similarity measures for indus-
reach their maximum level and the water level in tank 3 in- trial alarm flood analysis. In Industrial Informatics
creases. (INDIN), 2017 IEEE 15th International Conference
on, pages 781–786. IEEE, 2017.
[3] Raymond Reiter. A theory of diagnosis from first prin-
ciples. Artificial intelligence, 32(1):57–95, 1987.
[4] Elke Laubwald. Coupled tanks systems 1. control-
systems-principles. co. uk, 2015.
[5] Peter Struss. Fundamentals of model-based diagno-
sis of dynamic systems. In IJCAI (1), pages 480–485,
1997.
[6] Feng Lin. Diagnosability of discrete event systems
and its applications. Discrete Event Dynamic Systems,
4(2):197–212, 1994.
[7] Matthew J Daigle, Indranil Roychoudhury, Gautam
Biswas, Xenofon D Koutsoukos, Ann Patterson-Hine,
and Scott Poll. A comprehensive diagnosis method-
ology for complex hybrid systems: A case study on
spacecraft power distribution systems. IEEE Transac-
Figure 2: Development of the water level in tank 3 over tions on Systems, Man, and Cybernetics-Part A: Sys-
time in normal conditions (grey) and with a fault in valve 5 tems and Humans, 40(5):917–931, 2010.
(black) [8] Alban Grastien, Patrik Haslum, Sylvie Thiébaux, et al.
Conflict-based diagnosis of discrete event systems:
Theory and practice. In KR, 2012.
8 Conclusion [9] Indranil Roychoudhury, Matthew J Daigle, Gautam
Biswas, and Xenofon Koutsoukos. Efficient simula-
This paper shows for the limited domain of a fully-
tion of hybrid systems: A hybrid bond graph approach.
observable, hybrid, dynamic industrial production system
Simulation, 87(6):467–498, 2011.
that we can model its behaviour with state-space equations
and then translate it into satisfiability modulo theory and [10] Sriram Narasimhan and Gautam Biswas. Model-based
perform diagnosis. So far, operators in the process industry diagnosis of hybrid systems. IEEE Transactions on
rely only on fault identification techniques such as support- systems, man, and cybernetics-Part A: Systems and hu-
vector machines, artificial neural networks, Bayesian ap- mans, 37(3):348–361, 2007.
proaches etc. In this work we present an approach which can [11] Om Prakash and AK Samantaray. Model-based diag-
be used to also bring fault isolation into cyber-physical pro- nosis and prognosis of hybrid dynamical systems with
duction systems. Given a suitable model of the production dynamically updated parameters. In Bond Graphs for
system the presented method is able to capture behaviour Modelling, Control and Fault Diagnosis of Engineer-
over time while preserving the ability to directly react to ing Systems, pages 195–232. Springer, 2017.
faults. With a suitably chosen data sampling frequency the
[12] Alban Grastien. Diagnosis of hybrid systems by con-
potentially huge knowledge-base proposed by Grastien [8]
can be avoided. sistency testing. In 24th International Workshop on
For future work we will show how to extend this approach to Principles of Diagnosis (DX-13), pages 9–14. Cite-
deal with semi- and non-observable systems. For these we seer, 2013.
need better models of single components by, for example, [13] Martin Fränzle, Holger Hermanns, and Tino Teige.
creating assumables with the help of differential equations Stochastic satisfiability modulo theory: A novel tech-
specified in SMT logic. nique for the analysis of probabilistic hybrid systems.
Another direction for further research is the automatic gen- In International Workshop on Hybrid Systems: Com-
eration and learning of system models. Nowadays, models putation and Control, pages 172–186. Springer, 2008.
of systems need to be created manually which is not gener- [14] Hamed Khorasgani and Gautam Biswas. Structural
alizable and time consuming. The state-space equations and fault detection and isolation in hybrid systems. IEEE
their translations into logic are simple and can be specified Transactions on Automation Science and Engineering,
algorithmically. Therefore, an attempt should be made to 2017.
at least semi-automatically learn parts of these models from
descriptions of meta-data.
References
[1] Rolf Isermann and Peter Balle. Trends in the appli-
cation of model-based fault detection and diagnosis
of technical processes. Control engineering practice,
5(5):709–719, 1997.