# Security of differential phase shifted QKD against explicit individual
attacks

Valliamai Ramanathan, Anil Prabhakar, and Prabha Mandayam  Indian Institute of
Technology Madras, Chennai - 600036

(May 2023)

###### Abstract

Quantum key distribution (QKD) is known to be unconditionally secure in
principle, but quantifying the security of QKD protocols from a practical
standpoint continues to remain an important challenge. Here, we focus on
phase-based QKD protocols and characterize the security of the $3$ and
$n$-pulse Differential-Phase-Shifted Quantum Key Distribution (DPS QKD)
protocols against individual attacks. In particular, we focus on the minimum
error discrimination (MED) and cloning attacks and obtain the corresponding
bit error rates and the collision probability in the presence of these
attacks. We compare the secure key rates thus obtained with the known
theoretical lower bounds derived considering a general individual attack. In a
departure from the theoretical lower bounds which has no explicit attack
strategies, our work provides a practical assessment of the security of these
phase-based protocols based on attacks with known implementations.

###### Contents

  1. I Introduction
  2. II  Preliminaries
    1. II.1 Secure key rate
  3. III Minimum error discrimination attack
    1. III.1 MED of the ideal 3-pulse DPS QKD states
    2. III.2 Intercept and resend attack with MED
  4. IV Cloning attack
    1. IV.1 Cloning as an SDP problem
    2. IV.2 Optimal QCM for 3-pulse DPS states
    3. IV.3 Error introduced due to cloning
    4. IV.4 Collision probability and key rate with cloning attack
  5. V Conclusion

##  I Introduction

Quantum Key Distribution (QKD) offers the promise of secure communication over
public networks, ideally, with _unconditional security_ bb84 ; Bennett1992 ;
ekert1991quantum . In practice, this notion of unconditional security is
compromised by device imperfections and detector efficiencies. Furthermore,
different QKD implementations are susceptible to different kinds of
eavesdropping attacks, ranging from individual attacks like intercept-resend
IR_attack05  to photon number splitting (PNS) attacks pnsattack_95 ;
lutkenhaus2002quantum  and collective attacks biham2002security . We refer to
review articles practicalqkd2009 ; secureQKD2020  for a comprehensive survey
of the known secure key rate estimates, both in the ideal case as well as in
the imperfect scenarios.

The security of QKD relies on fundamental properties of quantum systems
including monogamy of entanglement, the no-cloning principle, and the
uncertainty principle shorpreskill . Quantitatively, the security of a
specific QKD implementation is characterised by its secure key rate Renner2008
; portmann2021security . Both asymptotic and finite key security analyses have
been presented for the standard QKD protocols in the literature
scarani2008quantum ; HayashiFinite ; sheridan2010finite ; cai2009finite . The
former is easier and assumes an infinite key length, which might be
unrealistic, whereas the latter takes into account the finite length size of
the key.

In this paper, we focus on the specific class of differential-phase-shifted
(DPS) QKD protocols Inoue ; InterceptResend2017 ; DPSInoue . This class of
phase-based protocols is known for its simplicity and efficient key
generation. DPS-QKD does not require basis reconciliation as in the case of
the BB84 protocol and thus every bit that is detected contributes to the key.
The ease of implementation of this protocol shaw2020equivalence  makes it
interesting from the implementation point of view. Several variations of this
protocol including round-robin DPS Sasaki2014RRDPS , the small-number-random
DPS protocol SNRDPS2017 , and measurement-device-independent DPS (MDI-DPS)
Ranu_2021  protocols have been studied in the recent past.

While the DPS protocol was originally proven to be secure against the basic
individual attacks such as the intercept-resend and beam splitter attacks, it
was subsequently shown to be secure against more general individual attacks
IndividualAttackSecurity . Unconditional security was then proved under
general coherent attacks in UnconditionalSecurity_PRL  for single photon DPS
QKD. In proving unconditional security, the eavesdropper is considered to be
capable of entangling ancillary systems to blocks of pulses and performing
operations that are dependent on previous measurement results. Other attacks
studied prior to the unconditional security bound UnconditionalSecurity_PRL ,
include specific collective attacks collecticeBS .

In this work, we study two specific kinds of individual attacks namely, the
minimum error discrimination (MED) attack and the cloning attack. The effects
of these attacks on the secure key rate of the DPS-QKD protocol still need to
be quantified. Note that the individual eavesdropping strategy, the simplest
class of strategies, is still technologically very challenging to achieve
today in an optimal way. The MED attack is similar to the unambiguous state
discrimination (USD) attack, except that the eavesdropper makes some errors in
identification. Recall that unambiguous state discrimination (USD) is only
possible for states that are linearly independent CHEFLES1998339 . Since the
ideal (single-photon) $n\geq 3$-pulse DPS states are linearly dependent, the
corresponding protocols are naturally secure against a USD-type attack. In
past work, InterceptResend2017 , an intercept and resend attack has been
considered wherein the eavesdropper (Eve) uses the same setup as the receiver
(Bob). In our work, we consider the scenario where Eve does minimum error
discrimination (MED) of the states after intercepting them. This constitutes
the first part of the study. MED is based on the detection theory devolved in
Helstrom1969  and has been studied in the past in the context of other phase-
based protocols such as the coherent one-way protocol MinErrorAttack .

The second part of our work analyzes the security of the DPS-QKD protocol
against a cloning attack. The cloning attack has been previously studied for
QKD protocols such as BB84 and the six-state protocol CloningOnBB84 ;
PhaseCovarintClonnerBB84 ; CloningAPS2004 , and has also been used as a
counterfeiting procedure for quantum money molina2012optimal . The problem of
finding the _optimal_ quantum cloner is a convex optimization problem and can
be formulated as a semidefinite program. In our work, we assume that the
eavesdropper uses a symmetric quantum cloning machine that maximizes the
cloning fidelity for the $3$-pulse DPS-QKD protocol. Both the MED measurement
and the cloning machine maps are optimized to achieve the maximum information
gain for an eavesdropper by recasting them as semidefinite programs.

The rest of the paper is organized as follows. In Sec. II we briefly review
the DPS-QKD protocol. Sec. III describes the MED of quantum states and
formulates the MED attack for the $3$-pulse and $n$-pulse DPS protocols. In
Sec. IV the description of a cloning attack against $3$-pulse DPS is provided,
along with the error rates introduced.

##  II  Preliminaries

We begin with a brief review of the $n$-pulse DPS-QKD protocol Inoue ;
InterceptResend2017 , schematically shown for the $n=3$ case in Fig. 1. We
note that this schematic can be easily extended to $n$ pulses by simply
introducing more delay lines at the source, or using an equivalent time-bin
phase modulation shaw2020equivalence . The sender (Alice) prepares a photon in
an equal superposition of $n\geq 3$ pulses and encodes her bit values
$\\{0,1\\}$ via the relative phases $\\{0,\pi\\}$. The states sent by Alice in
the $n=3$ case can be written as,

$|\psi\rangle_{\rm\pm,\pm}=\dfrac{1}{\sqrt{3}}(|100\rangle\pm|010\rangle\pm|001\rangle).$
(1)

Figure 1: The $3$-pulse DPS-QKD protocol. PM - Phase Modulator, $\theta_{12}$
and $\theta_{23}$ are the phase introduced between the pulses by the phase
modulator. ${\theta_{12},\theta_{23}\in\\{0,\pi\\}}$. D0 and D1 are single
photon detectors. Detection in the second and the third time slot gives
$\theta_{12}$ and $\theta_{23}$ respectively.

The receiver (Bob) uses an unbalanced Mach Zehnder interferometer (MZI) for
decoding, detects the relative phase between the pulses, and hence obtains the
bit value sent by Alice. On Bob’s side, the detection at the first and the
last time slots are random and are not used in generating the shared key.
Detection in any of the other $n-1$ time slots gives the key bits, wherein the
bit-value is assigned based on which detector (constructive or destructive)
clicks. In the ideal case, when a single photon is sent in a superposition of
$n$ pulses, the sifted key rate is obtained as,

$R_{\text{sifted}}=\frac{n-1}{n}.$ (2)

The security of this differential phase-shifted protocol is based on the fact
that the states in Eq. (1) are non-orthogonal and hence can not be
distinguished perfectly.

###  II.1 Secure key rate

The security of any QKD implementation can be quantified via the secure key
rate. Assuming that the eavesdropper is restricted to only individual attacks,
the secure key rate can be calculated for the DPS-QKD protocol, as DPSSecurity
,

$R_{\text{sk}}=R_{\rm{sifted}}\left\\{\tau-f(e)\left[h(e)\right]\right\\}.$
(3)

Here, $R_{\text{sifted}}=s\gamma p_{\text{click }}$, where $s$ is the sifting
parameter, with $s=2/3$ for $3$-pulse case. $\gamma$ is the repetition rate of
the transmission, $p_{\text{click }}$ is the probability of Bob’s detection
taking into account detector efficiencies, $e$ is the bit error rate, and
$f(e)$ characterizes the performance of the error correction algorithm. The
parameter $\tau$ represents the shrinking factor due to privacy amplification
and is calculated from the average collision probability and
$h(e)=-(e\log_{2}e+(1-e)\log_{2}(1-e))$ is the Shannon binary entropy
DPSSecurity . The key rate given above is in the asymptotic limit. The two
terms in the expression, $\tau$ and $f(e)[h(e)]$ relate to privacy
amplification and error correction respectively. Based on Shannon’s noiseless
coding theorem, the minimum number of bits that must be exchanged publicly for
error correction equals the binary entropy of the bit error rate $e$. The bit
error rate $e$ has contributions from both the signal and dark counts. The
error in the signal states can be due to imperfect state preparation. The dark
counts are uncorrelated to Alice’s key and thus have a 50 percent chance of
being wrong. The bit error rate is thus given as

$e=\dfrac{0.5p_{\text{dark}}+bp_{\text{signal}}}{p_{\text{click}}},$ (4)

where $b$ is the baseline system error rate quantifying the imperfect state
preparation. All the system’s errors are usually given as an advantage to Eve
and thus the role of privacy amplification is to deduce the shrinking factor
$\tau$ by which the key rate has to be shrunk in order to bound the
information leaked to Eve. The shrinking factor $\tau$ is a function of the
average collision probability $p_{c}$.

${\tau=-\log_{2}p_{c}}$ (5)

The average collision probability quantifies Eve’s mutual information with
Alice and Bob and can be evaluated in terms of the individual collision
probability of each bit, denoted as $p_{\text{co}}$. In fact, for individual
attacks, the overall collision probability is simply the product of the
individual collision probabilities. Specifically, $p_{\text{co}}$ is set to
$1$ for the bits that are detected correctly by Eve and $1/2$ for the bits
unknown to Eve DPSSecurity . In general, the average collision probability is
given asIndividualAttackSecurity

$p_{\text{co}}=\sum_{x,z}p^{2}(X_{i}=x|Z_{i}=z)p(Z_{i}=z)$ (6)

where X and Z are Alice’s and Eve’s key bit string and Eve’s information
respectively. In IndividualAttackSecurity , the most general individual attack
is considered and a bound on the collision probability for each bit is found
to be

$p_{\text{co}}\leq 1-e^{2}-{(1-6e)^{2}}/2.$

Using this value in Eq. (3) gives a lower bound on the secure key rate
obtained against individual attacks for DPS QKD using a weak coherent source
of mean photon number $\mu$ as IndividualAttackSecurity ,

${R_{\rm ind}=-p_{\rm click}[-(1-2\mu)\log_{2}p_{\rm co}+f(e)h(e)]}$ (7)

The term $(1-2\mu)$ is the fraction of bits that is completely known to Eve
due to a photon number splitting attack. The unconditional security proof in
UnconditionalSecurity_PRL  gives a bound on the phase error rate by assuming a
general coherent attack and hence a lower bound on the secure key rate is
derived to be

$R_{\text{sk}}\geq R_{\text{sifted}}[1-h(e)-h((3+\sqrt{5})e)].$ (8)

This expression holds for $n$-pulse DPS with $R_{\text{sifted}}$ depending on
$n$. In this paper, the adversary is said to apply an individual and identical
attack based on MED and cloning. There is no assumption that Eve holds any
information intercepted in quantum memory. We now proceed to describe the
attacks and the corresponding secure key rates.

##  III Minimum error discrimination attack

In the minimum error discrimination (MED) attack, the eavesdropper (Eve) aims
to identify the optimal quantum measurements required to discriminate the
signals transferred through a channel during QKD. The minimum error condition
is enforced by maximizing the success probability of the state discrimination.
MED is similar to the unambiguous state discrimination (USD) attack, except
that the eavesdropper is allowed to make some errors in identification. Recall
that unambiguous state discrimination (USD) is only possible for states that
are linearly independent CHEFLES1998339 . Since the ideal (single-photon)
$n$-pulse DPS states are linearly dependent for $n\geq 3$, the corresponding
protocols are naturally secure against a USD-type attack. Minimum error
discrimination can be described mathematically as follows - Suppose Alice has
a set of states $\\{\psi_{1},\psi_{2},...,\psi_{n}\\}$ which are non-
orthogonal and linearly dependent. She chooses $\psi_{i}$ with probability
$p(i)$ and sends it to Bob. Eve uses a set of positive operator valued
measures POVMs denoted as $\\{P_{1},...,P_{n}\\}$ to identify the states. The
POVM is set up such that $P_{i}$ when acting on the corresponding state
$\psi_{i}$ gives a result with high probability. In general, the average
success probability $\sum_{i}p(i)P_{i}\psi_{i}$ is maximized. Minimum error
discrimination can be stated as a semidefinite program (SDP), as described in
SDPResource . For completeness, we describe the SDP corresponding to the MED
as given in SDPResource . Let $\mathcal{X}$ be a complex Euclidean space,
$\Gamma$ a finite non-empty set, $p\in\mathbb{R}^{\Gamma}$ the probability
vector and $[\rho_{a}:a\in\Gamma]\subset D(\mathcal{X})$ a collection of
density operators. Alice randomly selects $a\in\Gamma$ according to $p$ and
sends the corresponding $\rho_{a}$ to Eve. Eve uses the measurement $\mu$,
$[\mu:\Gamma\rightarrow\text{Pos}(\mathcal{X}):a\rightarrow P_{a}]$.

The SDP for performing MED is given as

$\begin{array}[]{ll}\text{maximize}&\sum_{a\in\Gamma}\langle\sigma_{a}P_{a}\rangle\\\
\text{ subject to
}&\\{P_{a}:a\in\Gamma\\}\subset\operatorname{Pos}(\alpha),\\\
&\sum_{a\in\Gamma}P_{a}=\mathds{1}_{\mathcal{X}}\end{array}$ (9)

where $\sigma_{a}=p(a)\rho_{a}$, the product of the density matrix sent by
Alice and the probability with which she sends it.$P_{a}$ are the POVM
elements used by Eve to perform the attack. Thus a SDP solver can be used to
find the success probabilities and the POVM operators for minimum error
discrimination of quantum states.

###  III.1 MED of the ideal 3-pulse DPS QKD states

The three-pulse ideal single photon DPS QKD shown in Fig.1 encodes the logical
bits in the phase information. The states sent by Alice can be written as in
Eq. (1). We assume that Alice sends the states with equal probability i.e.
$\frac{1}{4}$. Thus to do a minimum error discrimination of these four states,
Eve has to set up a POVM set of four elements -
$\\{P_{1},P_{2},P_{3},P_{4}\\}$. Each POVM element rightly identifies one
state. The states in vector form can be written as,

$\begin{array}[]{ll}\psi(+,+)=(1/\sqrt{3})[0;1;1;0;1;0;0;0];\\\
\psi(+,-)=(1/\sqrt{3})[0;-1;1;0;1;0;0;0];\\\
\psi(-,+)=(1/\sqrt{3})[0;1;-1;0;1;0;0;0];\\\
\psi(-,-)=(1/\sqrt{3})[0;-1;-1;0;1;0;0;0].\end{array}$

The SDP corresponding to the DPS states minimum error discrimination is-

$\begin{array}[]{ll}\text{
maximize:}&\frac{1}{4}\langle\rho_{i},P_{i}\rangle\\\ \text{subject
to:}&\sum_{i}P_{i}=I\\\ &P_{i}>0\end{array}$ (10)

Here, $\rho_{i}$ is the density matrix corresponding to the DPS states.

The SDP is run using the cvx solver cvx ; gb08  and the result is summarised
in Table 1. We see that each state has a 75% chance of being rightly
identified by Eve.

Table 1: Measurement probabilities for the $3$-pulse DPS states | $P_{1}$ | $P_{2}$ | $P_{3}$ | $P_{4}$
---|---|---|---|---
$\psi(++)$ | $0.75$ | $0.083$ | 0.083 | 0.083
$\psi(+-)$ | 0.083 | $0.75$ | 0.083 | 0.083
$\psi(-+)$ | 0.083 | 0.083 | 0.75 | 0.083
$\psi(--)$ | 0.083 | 0.083 | 0.083 | 0.75

###  III.2 Intercept and resend attack with MED

The MED of the 3-pulse DPS states shows that Eve has a 0.75 probability of
rightly identifying a given state. In InterceptResend2017  an intercept and
resend attack was considered where Eve uses the same apparatus as Bob, and Eve
introduces errors in 33% of the intercepted bits. In an intercepting attack
using MED, Eve introduces only a 25% error. We now compare the key rate in the
presence of these two attacks. In order to find the secure key rate in the
presence of the minimum error discrimination attack, one has to find the
collision probability and the corresponding shrinking factor. The collision
probability is defined as in eq. (6), with $P(X|Z)$ being the probability of
guessing Alice’s string given Eve’s information Z. In the case of the MED
attack, Eve’s information corresponds to which of the POVMs clicked. Therefore
$Z\in\\{P_{1},P_{2},P_{3},P_{4}\\}$. Alice’s bit value is either 0 or 1 i.e
$X\in\\{0,1\\}$, which can in turn be detected from the phase difference
between the first and the second time bin or the second and third, both the
cases being equally probable. Now,

$\displaystyle p_{\rm{co}}=$ $\displaystyle\sum_{x,z}P^{2}(X=x|Z=z)P(Z=z)$
$\displaystyle=$ $\displaystyle\sum_{i=1}^{4}P^{2}(X=0|Z=P_{i})P(P_{i})+$
$\displaystyle P^{2}(X=1|Z=P_{i})P(P_{i})$ $\displaystyle=$ $\displaystyle
1/2[~{}\sum_{i=1}^{4}P^{2}(X=0_{1}|Z=P_{i})P(P_{i})+$ (11) $\displaystyle
P^{2}(X=1_{1}|Z=P_{i})P(P_{i})~{}]$
$\displaystyle+1/2[~{}\sum_{i=1}^{4}P^{2}(X=0_{2}|Z=P_{i})P(P_{i})+$
$\displaystyle P^{2}(X=1_{2}|Z=P_{i})P(P_{i})]$

The final step is to account for the fact that the bit detected can correspond
to either the first phase difference or the second. The probabilities required
to calculate eq(III.2) are given in table 1 and the average collision
probability is found to be 0.72. Now, Let the probability of correct detection
by Eve be denoted as $p_{\text{corr}}$, and the fraction of bits intercepted
by Eve as $n_{\text{int}}^{\text{Eve}}$. To maintain the error rate at Bob’s
side, the fraction of bits Eve intercepts is,

$(1-p_{\text{corr}})n_{\text{int}}^{\text{Eve}}=e.$ (12)

Now, from the SDP, we have found that $p_{\text{corr}}=0.75$ for the minimum
error discrimination of the $3$-pulse DPS states. This indicates that the
fraction of bits that can be intercepted by Eve is $4e$. Out of these 4e
pulses intercepted only two-thirds of 4e is detected by bob in either the
second or in the third time bin. Thus for a fraction of $(2/3)*4e$ the
collision probability is 0.72 while for the rest ($1-(2/3)4e$) the collision
probability is half since it is unknown to Eve. The key rate in the presence
of MED is thus given as

$R_{\rm MED}=-s\gamma p_{\rm click}[-((2/3)4e)\log_{2}0.72\\\
+(1-(2/3)4e)+f(e)h(e)]$ (13)

The key rate is calculated for the following implementation parameters:
$\gamma=10$ MHz, channel loss of 0.2 dB/km, $f(e)$ 1.16. The key rate is
plotted in Fig. 2 and it can be seen that the length up to which the protocol
is secure reduces in the presence of the intercept and resend attack with MED
as compared to that with intercept and resend using an MZI. The SDP result and
key rate plots can be recreated using the codes in plot_git .

Figure 2: Secure Key rate vs channel length

In the four pulse case, the number of different combinations of the phases is
8. Thus here Eve has to discriminate between 8 different states to rightly
resend.

The states are

$|\psi\rangle_{\rm\pm,\pm,\pm}=\dfrac{1}{\sqrt{4}}(|1000\rangle\pm|0100\rangle\pm|0010\rangle\pm|0001\rangle)$
(14)

Now, a POVM of 8 elements is set up, and the probability of the right
identification is maximized. This SDP problem gives out that each state has a
probability of 0.5 to be rightly identified. That is in every attempt she just
has a 50% chance of right identification. Thus in the four-pulse case, the
minimum error attack is not a strong attack. Similarly, in the 5 pulse case,
the number of states goes to 16 and the probability of right identification
drops further down to 0.31. Thus identification becomes erroneous 75% of the
time. From figure 3 it can be noted that the minimum error discrimination
attack using this formulation is not powerful in the case of $n>3$ pulses.

Figure 3: Fraction of bits erroneous due to MED as $n$ increases

##  IV Cloning attack

In this attack, the eavesdropper has a Quantum Cloning Machine (QCM) and tries
to clone the photon state sent to Bob. Since the theory of quantum mechanics
does not allow for perfect cloning Wootters1982 , the cloned state is not a
perfect copy and the initial state is disturbed to some level. In this
section, we find the bit error rate due to such imperfect cloning.

###  IV.1 Cloning as an SDP problem

The cloning attack is described by a quantum channel $\Phi$, which takes a
state $\rho=|\psi\rangle\langle\psi|\in\mathrm{D}(\mathcal{X})$ to the state
$\Phi(\rho)\in\mathrm{D}(\mathcal{Y}\otimes\mathcal{Z})$ where $\mathrm{D}$
refers to the set of density matrices in the Hilbert space specified. The
channel $\Phi$ must correspond to a completely positive and trace-preserving
linear mapping of the form
$\Phi:\mathrm{L}(\mathcal{X})\rightarrow\mathrm{L}(\mathcal{Y}\otimes\mathcal{Z}).$
Suppose the state that is sent is indexed $k$, the probability of success for
an attack described by $\Phi$ is given by
$\left\langle\psi_{k}\otimes\psi_{k}\left|\Phi\left(\left|\psi_{k}\right\rangle\left\langle\psi_{k}\right|\right)\right|\psi_{k}\otimes\psi_{k}\right\rangle.$
Averaging over the possible choices of $k$, the overall success probability of
a cloning attack is

$\sum_{k=1}^{N}p_{k}\left\langle\psi_{k}\otimes\psi_{k}\left|\Phi\left(\left|\psi_{k}\right\rangle\left\langle\psi_{k}\right|\right)\right|\psi_{k}\otimes\psi_{k}\right\rangle$
(15)

The optimal success probability of a cloning strategy, which is represented by
the supremum of the probability over all valid channels
$\Phi:\mathrm{L}(\mathcal{X})\rightarrow\mathrm{L}(\mathcal{Y}\otimes\mathcal{Z})$,
can be represented by a semidefinite program.

The formulation makes use of the Choi-Jamiołkowski representation
$J(\Phi)\equiv X$ of a given channel $\Phi$. The SDP is given as
molina2012optimal

$\begin{array}[]{ll}\text{maximize}&\langle Q,X\rangle\\\ \text{ subject
to}&\text{Tr}_{\mathcal{Y}\otimes\mathcal{Z}}(X)=\mathbb{I}_{\mathcal{X}}\\\
&\mathcal{X}\in\operatorname{Pos}(\mathcal{Y}\otimes\mathcal{X}\otimes\mathcal{Z})\end{array}$
(16)

where

$Q=\sum_{k=1}^{N}p_{k}\left|\psi_{k}\otimes\psi_{k}\otimes\overline{\psi_{k}}\right\rangle\left\langle\psi_{k}\otimes\psi_{k}\otimes\overline{\psi_{k}}\right|$

Here, $\overline{\psi}$ refers to taking complex conjugation with respect to
the standard basis.

###  IV.2 Optimal QCM for 3-pulse DPS states

To find the cloning map for the DPS states, let us consider the DPS three
pulse state to be the states of a qutrit with basis states
$\\{|1\rangle,|2\rangle,|3\rangle\\}$. The four DPS states can be then written
as,

$|\psi\rangle_{\rm\pm,\pm}=\dfrac{1}{\sqrt{3}}(|1\rangle\pm|2\rangle\pm|3\rangle)$
(17)

A cloning machine to clone this set of qutrit states can be divised using Eq
(16). The SDP solution for the probability of right cloning of the qutrit
state gives a value of 0.78 plot_git . We can obtain the map acting on Bob’s
state from the Choi matrix. The fidelity $\mathcal{F}$ between the original
state and the cloned state is 0.81. This fidelity is obtained for all the
states $|\psi(\pm,\pm)\rangle$.

The bit error rate at Bob’s side due to the application of a cloning machine
is found out using the probability of detecting in the constructive
(destructive ) port when the destructive (constructive) port is supposed to
click. We write Bob’s mixed state density matrix after the cloning machine as
a convex combination on an appropriate basis and then find out the probability
of wrong detection.

###  IV.3 Error introduced due to cloning

In the three pulse DPS protocol, the states sent by Alice are

$|\psi\rangle_{\rm\pm,\pm}=\dfrac{1}{\sqrt{3}}(|1\rangle\pm|2\rangle\pm|3\rangle)$
(18)

This is after writing the DPS 3-pulse state as a qutrit state, the basis
states here are $|1\rangle$ which is equivalent to the photon in the third
time bin $|100\rangle$, $|2\rangle$ which is equivalent to the photon in the
second time bin $|010\rangle$, $|3\rangle$ which is equivalent to the photon
in the first time bin $|001\rangle$. The input states and the Choi matrix for
the SDP are written in the qutrit basis. The SDP gives the Choi matrix from
which the map acting on Bob’s state is obtained. For the input,

$\rho(++)=\frac{1}{3}\left(\begin{array}[]{ccc}1&1&1\\\ 1&1&1\\\
1&1&1\end{array}\right)$ (19)

After cloning, the reduced density matrix reaching Bob and Eve is,

$\left(\begin{array}[]{ccc}\dfrac{1}{3}&0.23&0.23\\\ 0.23&\dfrac{1}{3}&0.23\\\
0.23&0.23&\dfrac{1}{3}\end{array}\right)$ (20)

An error occurs in Bob’s detector if there is a detection in the wrong port.
Let us take the example of the state in (19). The detector at Bob’s side is a
1-bit delay MZI with the input-output relation as shown in Fig. 4.

Figure 4: 1-bit delay MZI

Now, when (19) is sent, the state can be written as

$\psi(++)=\frac{1}{\sqrt{3}}\left[a_{\rm k+2}^{\dagger}+a_{\rm
k+1}^{\dagger}+a_{\rm k}^{\dagger}\right]|0\rangle$ (21)

After $\psi(++)$ goes through the MZI, the output state is,

$\frac{1}{\sqrt{3}}\left[\frac{u_{\rm k}^{\dagger}}{2}+\frac{iv_{\rm
k}^{\dagger}}{2}+u_{\rm k+1}^{\dagger}+u_{\rm k+2}^{\dagger}+\frac{u_{\rm
k+3}^{\dagger}}{2}-\frac{iv_{\rm k+3}^{\dagger}}{2}\right]|0\rangle$ (22)

In the $k+1$ and $k+2$ slots, only the constructive port clicks. Now let us
look at what happens when the state sent by Alice is cloned and the state
received by Bob is given by (20). To look at the probabilities of detection in
different ports, we write the mixed state in terms of a basis state given as
$[\frac{1}{\sqrt{3}}(|1\rangle+|2\rangle+|3\rangle),\frac{1}{\sqrt{6}}\left(|1\rangle+|2\rangle-2|3\rangle\right),\frac{1}{\sqrt{2}}\left(|1\rangle-|2\rangle\right)]$.
Denoting the basis density matrices as $\rho_{1}$, $\rho_{2}$, $\rho_{3}$. The
mixed-density matrix received by Bob can be written as,

$0.81\rho_{1}+0.095\rho_{2}+0.095\rho_{3}$ (23)

The second and third density operators can easily be seen to cause a detection
in the destructive port. This erroneous click is observed due to $\rho_{2}$ at
$k+2$ slot and $\rho_{3}$ in both slots $k+1$ and $k+2$. After passing the
cloned state through the MZI at Bob’s side the probability of error (in this
case - the probability of clicking in the destructive port in the second and
third-time bin) is calculated as:

$\text{BER}=0.095\left(\frac{3}{8}\right)+0.095=0.13$ (24)

Each cloned state has a probability of giving a wrong detection of 0.13. Thus
an eavesdropper in order to hide his/her presence will try to clone only some
number of states such that the error introduced matches the error introduced
due to the channel.

$p_{\text{error}}n_{\text{int}}^{\text{Eve}}=e$

If the error rate in the system from signal and dark counts is $e$, for the
error probability of $13\%$ Eve can intercept $7.69e$ fraction of bits and
remain unidentified.

###  IV.4 Collision probability and key rate with cloning attack

In order to quantify the maximum information available to an adversary, we let
Eve do a minimum error discrimination of the states reaching her after the
cloning machine. Since the states reaching Eve and Bob are mixed, the MED
success probability is relatively less than in the case where Eve does a MED
directly on the states sent by Alice. Nevertheless, this strategy does
marginally better since the error introduced at Bob due to cloning is less
than the direct MED thus allowing for more bits to be attacked. Denoting the
POVM elements for MED of the cloned state as $\rm{CP_{i}}$, the SDP results
for the right identification are as given in table 2. The collision
probability is found for these POVM as in the case of the MED attack and is
found to be 0.61. The corresponding secure key rate is given as -

$R_{\rm Cloning}=-s\gamma p_{\rm click}[-((2/3)7.69e)\log_{2}0.61+\\\
(1-(2/3)7.69e)+f(e)h(e)]$ (25)

Figure 5: Secure Key rate vs channel length  Table 2: Measurement probabilities for the $3$-pulse DPS states after cloning | $\rm{CP_{1}}$ | $\rm{CP_{2}}$ | $\rm{CP_{3}}$ | $\rm{CP_{4}}$
---|---|---|---|---
$\psi(++)$ | $0.6071$ | $0.1310$ | 0.1310 | 0.1310
$\psi(+-)$ | 0.1310 | $0.6071$ | 0.1310 | 0.1310
$\psi(-+)$ | 0.1310 | 0.13103 | 0.6071 | 0.1310
$\psi(--)$ | 0.1310 | 0.1310 | 0.1310 | 0.6071

Finally, we plot the secure key rate vs the distance for all the attacks
discussed so far in Fig. 5

Table 3: Bit Error rates introduced in various attacks. $\mu$ \- mean photon number, T - transmission loss, e - error rate due to implementation. The probability of success is the probability of the right identification of a state. I & R refer to intercept and resend attack Protocol | General I & R | I & R - USD/MED | Cloning attack
---|---|---|---
Pulse train DPSInoue  (WCS) |  BER = 0.25, 2e fraction with Eve DPSSecurity  |  USD prob. of success $=1-e^{-2\mu}$ | -
3-pulse DPS (Single photon) |  BER = 0.33 e fraction with Eve InterceptResend2017  |  MED prob. of success = 0.75 |  prob. of success = 0.78, BER = 0.13
3-pulse DPS (WCS) |  BER= 0.33, $\frac{2}{9}e^{-\mu}\mu$ fraction with Eve |  USD prob. of success $=1-e^{-2\mu}$ | -
n-pulse DPS | Learning rate = $\frac{n-1}{n^{2}}$Ranu_2021  |  MED $n=4,p_{\text{succ}}=0.50$ $n=5,p_{\text{succ}}=0.31$ | -

##  V Conclusion

We quantify the secure key rate for the $3$-pulse DPS QKD protocol in the
presence of minimum error discrimination and cloning attack. For both these
attacks, we have identified the optimal strategy for the eavesdropper, for the
specific set of signal states used in the $3$-pulse DPS protocol. We further
note that these attack strategies can be implemented by realising specific
maps and are realistic in the sense that they do not require quantum memory.

Our results indicate that these sophisticated attacks are in fact comparable
to the simple intercept and resend attack in terms of the information gained
by the adversary. An eavesdropper gets only as much information as a simple
”Bob-like” intercept and resend strategy even if she uses the minimum error
discrimination or cloning attack as described here. These explicit attack
strategies give much higher secure key rates than the lower bound predicted
using inequalities such as the Cauchy inequality. As a point of comparison, we
summarize the various known individual attacks and the corresponding figures
of merit in Table 3.

Our work fills a certain gap that existed in the literature, with regards to
individual attacks on the DPS QKD protocol. However, we have entirely focused
on the ideal, single-photon DPS protocol here. Going forward, it would be
interesting to see how these individual attacks fare in the context of more
realistic phase-based QKD protocols involving weak coherent sources
coherentcloning .

Although unconditional security of DPS-QKD is known under general, collective
attacks, it is useful to benchmark the effectiveness of specific individual
attacks so as to be able to assess their relative strengths. The fact that the
secure key rates in the presence of the individual attacks studied here are
still much higher than the lower bound on the secure key rate also suggests
that these rates are of more practical relevance. Essentially, the secure key
rates under such explicit attacks give a better sense of the realistic key
rates and distances for QKD protocols.

## References

  * (1) C. H. Bennett and G. Brassard, “Quantum cryptography: Public key distribution and coin tossing,” in Proceedings of IEEE International Conference on Computers, Systems, and Signal Processing, p. 175â€“179, December 1984.
  * (2) C. H. Bennett, F. Bessette, G. Brassard, L. Salvail, and J. Smolin, “Experimental quantum cryptography,” Journal of Cryptology, vol. 5, pp. 3–28, Jan 1992.
  * (3) A. K. Ekert, “Quantum cryptography based on Bell’s theorem,” Phys. Rev. Lett., vol. 67, pp. 661–663, Aug 1991.
  * (4) M. Curty and N. Lütkenhaus, “Intercept-resend attacks in the Bennett-Brassard 1984 quantum-key-distribution protocol with weak coherent pulses,” Phys. Rev. A, vol. 71, p. 062301, Jun 2005.
  * (5) B. Huttner, N. Imoto, N. Gisin, and T. Mor, “Quantum cryptography with coherent states,” Phys. Rev. A, vol. 51, pp. 1863–1869, Mar 1995.
  * (6) N. Lütkenhaus and M. Jahma, “Quantum key distribution with realistic states: photon-number statistics in the photon-number splitting attack,” New Journal of Physics, vol. 4, no. 1, p. 44, 2002.
  * (7) E. Biham, M. Boyer, G. Brassard, J. Van De Graaf, and T. Mor, “Security of quantum key distribution against all collective attacks,” Algorithmica, vol. 34, no. 4, pp. 372–388, 2002.
  * (8) V. Scarani, H. Bechmann-Pasquinucci, N. J. Cerf, M. Dušek, N. Lütkenhaus, and M. Peev, “The security of practical quantum key distribution,” Rev. Mod. Phys., vol. 81, pp. 1301–1350, Sep 2009.
  * (9) F. Xu, X. Ma, Q. Zhang, H.-K. Lo, and J.-W. Pan, “Secure quantum key distribution with realistic devices,” Rev. Mod. Phys., vol. 92, p. 025002, May 2020.
  * (10) P. W. Shor and J. Preskill, “Simple proof of security of the BB84 quantum key distribution protocol,” Phys. Rev. Lett., vol. 85, pp. 441–444, Jul 2000\.
  * (11) R. RENNER, “Security of quantum key distribution,” International Journal of Quantum Information, vol. 06, no. 01, pp. 1–127, 2008.
  * (12) C. Portmann and R. Renner, “Security in quantum cryptography,” Rev. Mod. Phys., vol. 94, p. 025008, Jun 2022.
  * (13) V. Scarani and R. Renner, “Quantum cryptography with finite resources: Unconditional security bound for discrete-variable protocols with one-way postprocessing,” Physical review letters, vol. 100, no. 20, p. 200501, 2008\.
  * (14) M. Hayashi, “Upper bounds of eavesdropper’s performances in finite-length code with the decoy method,” Phys. Rev. A, vol. 76, p. 012329, Jul 2007.
  * (15) L. Sheridan, T. P. Le, and V. Scarani, “Finite-key security against coherent attacks in quantum key distribution,” New Journal of Physics, vol. 12, no. 12, p. 123019, 2010.
  * (16) R. Y. Cai and V. Scarani, “Finite-key analysis for practical implementations of quantum key distribution,” New Journal of Physics, vol. 11, no. 4, p. 045024, 2009.
  * (17) K. Inoue, E. Waks, and Y. Yamamoto, “Differential phase shift quantum key distribution,” Phys. Rev. Lett., vol. 89, p. 037902, Jun 2002.
  * (18) S. Ranu, G. Shaw, A. Prabhakar, and P. Mandayam, “Security with 3-pulse differential phase shift quantum key distribution,” 2017 IEEE Workshop on Recent Advances in Photonics (WRAP), 2017.
  * (19) K. Inoue, “Differential phase-shift quantum key distribution systems,” IEEE Journal of Selected Topics in Quantum Electronics, vol. 21, no. 3, pp. 109–115, 2015.
  * (20) G. Shaw, S. Sridharan, S. Ranu, F. Shingala, P. Mandayam, and A. Prabhakar, “Time-Bin Superposition Methods for DPS-QKD,”, IEEE Photonics Journal, vol. 14, no. 5, pp. 1-7 ,2022.
  * (21) T. Sasaki, Y. Yamamoto, and M. Koashi, “Practical quantum key distribution protocol without monitoring signal disturbance,” Nature, vol. 509, pp. 475–478, May 2014.
  * (22) Y. Hatakeyama, A. Mizutani, G. Kato, N. Imoto, and K. Tamaki, “Differential-phase-shift quantum-key-distribution protocol with a small number of random delays,” Phys. Rev. A, vol. 95, p. 042301, Apr 2017.
  * (23) S. K. Ranu, A. Prabhakar, and P. Mandayam, “Differential phase encoded measurement-device-independent quantum key distribution,” Quantum Information Processing, Volume 20, Article no. 67, Feb 2021.
  * (24) E. Waks, H. Takesue, and Y. Yamamoto, “Security of differential-phase-shift quantum key distribution against individual attacks,” Phys. Rev. A, vol. 73, p. 012344, Jan 2006.
  * (25) K. Wen, K. Tamaki, and Y. Yamamoto, “Unconditional security of single-photon differential phase shift quantum key distribution,” Phys. Rev. Lett., vol. 103, p. 170503, Oct 2009.
  * (26) C. Branciard, N. Gisin, and V. Scarani, “Upper bounds for the security of two distributed-phase reference protocols of quantum cryptography,” New Journal of Physics, vol. 10, p. 013031, jan 2008.
  * (27) A. Chefles, “Unambiguous discrimination between linearly independent quantum states,” Physics Letters A, vol. 239, no. 6, pp. 339–347, 1998.
  * (28) C. W. Helstrom, “Quantum detection and estimation theory,” Journal of Statistical Physics, vol. 1, pp. 231–252, Jun 1969.
  * (29) C. Branciard, N. Gisin, N. Lutkenhaus, and V. Scarani, “Zero-error attacks and detection statistics in the coherent one-way protocol for quantum cryptography,” Quantum Info. Comput., vol. 7, p. 639-664, sep 2007.
  * (30) L. Gyongyosi and S. Imre, “Computational geometric analysis of physically allowed quantum cloning transformations for quantum cryptography,” in Proceedings of the 4th WSEAS International Conference on Computer Engineering and Applications, CEA’10, (Stevens Point, Wisconsin, USA), p. 121-126, World Scientific and Engineering Academy and Society (WSEAS), 2010.
  * (31) Y. Okubo, F. Buscemi, and A. Tomita, “Proposal of an eavesdropping experiment for BB84 QKD protocol with 1-3 phase covariant quantum doner,” AIP Conference Proceedings, vol. 1110, no. 1, pp. 355–358, 2009.
  * (32) T. Durt, D. Kaszlikowski, J.-L. Chen, and L. C. Kwek, “Security of quantum key distributions with entangled qudits,” Phys. Rev. A, vol. 69, p. 032313, Mar 2004.
  * (33) A. Molina, T. Vidick, and J. Watrous, “Optimal counterfeiting attacks and generalizations for Wiesner’s quantum money,” 2012.
  * (34) E. Diamanti, “Security and implementation of differential phase shift quantum key distribution systems,” ProQuest Dissertations and Theses, p. 166, 2006\.
  * (35) J. Watrous, “Lecture notes on quantum information theory - lecture 1,2,7,8.”
  * (36) C. R. Inc.”, “Cvx: Matlab software for disciplined convex programming,” https://cvxr.com/cvx, 2011.
  * (37) M. Grant and S. Boyd, Graph implementations for nonsmooth convex programs.  Lecture Notes in Control and Information Sciences, Springer-Verlag Limited, 2008.
  * (38) GitHub, “To recreate the secure key rate plots.” https://github.com/ValliamaiRM/DPS_Keyrate, 2023.
  * (39) W. K. Wootters and W. H. Zurek, “A single quantum cannot be cloned,” Nature, vol. 299, pp. 802–803, Oct 1982.
  * (40) Andersen, Ulrik L. and Josse, Vincent and Leuchs, Gerd, “Unconditional Quantum Cloning of Coherent States with Linear Optics,” Phys. Rev. Lett., vol. 94, issue 24, 2005.


11institutetext: Department of Physics, Fudan University, Shanghai, China
22institutetext: Digital Medical Research Center, School of Basic Medical
Sciences, Fudan University, Shanghai, China 33institutetext: Shanghai Key
Laboratory of Medical Image Computing and Computer Assisted Intervention,
Shanghai, China 44institutetext: Shanghai Zhiyu Software Information Co., Ltd
55institutetext: Human Phenome Institute, Fudan University, Shanghai, China
66institutetext: Institute of Science and Technology for Brain-inspired
Intelligence, Fudan University, Shanghai, China 77institutetext: Key
Laboratory of Computational Neuroscience and Brain-Inspired Intelligence
(Fudan University), Ministry of Education, China
77email<EMAIL_ADDRESS>
http://homepage.fudan.edu.cn/hewang

# Non-invasive Assessment of Hepatic Venous Pressure Gradient (HVPG) Based on
MR Flow Imaging and Computational Fluid Dynamics

Kexin Wang 11 Shuo Wang 2233 Minghua Xiong 44
Chengyan Wang 55 He Wang 556677

###### Abstract

Clinically significant portal hypertension (CSPH) is a severe complication of
chronic liver disease associated with cirrhosis, which is diagnosed by the
measurement of hepatic venous pressure gradient (HVPG). However, HVPG
measurement is invasive and therefore difficult to be widely applied in
clinical routines. There is no currently available technique to measure HVPG
noninvasively. Computational fluid dynamics (CFD) has been used for
noninvasive measurement of vascular pressure gradient in the intracranial and
coronary arteries. However, it has been scarcely employed in the hepatic
vessel system due to the difficulties in reconstructing precise vascular
anatomies and setting appropriate boundary conditions. Several computer
tomography and ultrasound based studies have verified the effectiveness of
virtual HVPG (vHVPG) by directly connecting the portal veins and hepatic veins
before CFD simulations [16, 12]. We apply the latest techniques of phase-
contrast magnetic resonance imaging (PC-MRI) and DIXON to obtain the velocity
and vessel anatomies at the same time. Besides, we improve the CFD pipeline in
regards to the construction of vessel connections and reduction of calculation
time. The proposed method shows high accuracy in the CSPH diagnosis in a study
containing ten healthy volunteers and five patients. The MRI-based noninvasive
HVPG measurement is promising in the clinical application of CSPH diagnosis.

###### Keywords:

Computational fluid dynamics Hepatic venous pressure gradient Phase-contrast
MRI.

##  1 Introduction

Portal hypertension is a severe complication among patients with chronic liver
disease. Currently, clinically significant portal hypertension (CSPH) is
diagnosed with the measurement of hepatic venous pressure gradient (HVPG)
larger than 10 mmHg. This HVPG measurement is conducted by putting a balloon
catheter into the right hepatic vein to get both the free hepatic venous
pressure (FHVP) and the wedged hepatic venous pressure (WHVP). The HVPG is
defined as the difference between WHVP and FHVP. Although the gold standard of
clinical diagnosis of CSPH is HVPG, it is invasive, probably leading to
infection[6, 7] and other complications [4], thus inaccessible to some
patients.

Recently, several non-invasive measurements have been proposed to obtain HVPG,
transient elastography (TE) for example[14, 1]. However, TE is often
unreliable for patients with obesity or intrahepatic inflammatory activities
[3]. Some serum biomarkers including prothrombin index (PI) and aspartate
aminotransferase (AST) to alanine aminotransferase (ALT) ratio (AAR) are more
widely available in hospitals, but their accuracy is lower than TE.

Computational fluid dynamics (CFD) simulation has grown rapidly in the past
few years for the diagnosis of cardiovascular diseases[18, 13, 2], carotid
artery diseases[10] and intracranial artery occlusive diseases[9]. The virtual
blood pressure gradient is measured by solving the Navier-Stokes equation with
certain boundary conditions. CFD-based measurement of blood pressure gradient
has achieved remarkable success in assessing intracranial and coronary
arteries (e.g., cardiac fractional flow reserve (FFR)), while it is less
explored in the hepatic vascular system due to two possible reasons. First,
accurate CFD simulation is highly dependant on the precise reconstruction of
vascular anatomy and the measurement of inflow velocity. However, the
anatomies of hepatic vein and portal vein are difficult to be obtained.
Currently available ultrasound (US) based velocity measurement is not robust
and repeatable. Second, the computational efficiency is another limitation for
its application in clinical routines due to the complex structure of liver
vessels.

This paper proposes a one-stop solution for the noninvasive measurement of
HVPG in both the healthy volunteers and patients with liver cirrhosis by using
multi-contrast magnetic resonance imaging (MRI) and an efficient CFD model.

##  2 Methods and Experiments

###  2.1 Multi-contrast MRI

To avoid the use of contrast agents, multi-echo DIXON (mDIXON) imaging is
conducted to obtain the anatomy of hepatic vessels in this study. A 3D
gradient-echo (GRE) sequence with multiple echoes is used. All the subjects
are scanned on a 3.0 T MR scanner (Ingenia, Philips, the Netherlands). The
imaging parameters are: repetition time (TR) = 3.4 ms; echo time (TE) = 0.95
ms; field-of-view (FOV) = 350 × 350 × 120 mm3; flip angle (FA) = 10°; spatial
resolution = 0.87 × 0.87 × 2.0 mm3, which sets the limit for the smallest
vessel reconstructed; SENSitivity Encoding (SENSE) factor = 2; scan time =
12.5 ms.

Blood flow velocity is measured using phase-contrast (PC) GRE sequence.
Imaging is performed during free breathing with cardiac and respiratory
gatings. The imaging parameters for PC-MRI are: TR = 2.9 ms; TE = 1.45 ms; FOV
= 300 × 300 × 300 mm3; FA = 6°; spatial resolution = 0.87 × 0.87 × 2.0 mm3;
SENSE factor = 2; encoding velocity = 40 cm/s. The acquisition time is nearly
2.5 minutes depending on the subject’s heart rate. The blood flow acquisition
plane is selected perpendicular to the portal vein on the morphological scout
sequences.

###  2.2 Image Analysis

MRI data are transferred in Digital Imaging and Communications in Medicine
(DICOM) format to a workstation for analysis. The vessels are automatically
extracted from the mDIXON images using an pre-trained deep neural networks
(DNN)[17]. Two separate pre-trained deep neural networks (DNNs) are applied to
obtain the portal veins and hepatic veins automatically. The DNNs are trained
on an additional dataset including 30 manually labelled cases in our previous
work.

Extraction of the temporal flow velocity curve during a cardiac cycle is from
the cross section of portal vein manually labeled by an experienced
radiologist on the PC-MRI images. The flow velocity is calculated as the mean
value of the cross section obtained during the 30 phases of several repeated
cardiac cycles (see the process in Fig. 1. $\boldsymbol{u}(t)$).

Figure 1: The schematic diagram of the MRI based CFD model. From the MRI
images on the left to the pressure distribution at the bottom right. Firstly,
geometry of the hepatic vein and the portal vein are segmented individually,
together with the attained velocity of flow at the entrance of the portal
vein. The combination of the yellow and blue arrow is the green arrow, which
represents applying the connected model and boundary conditions, then with
meshing based on both practice, CFD calculation achieves the pressure
distribution eventually.

###  2.3 CFD Simulation

CFD simulation can provide accurate pressure measurement based on the geometry
of blood vessels and the boundary conditions, e.g., the velocity and the
reference pressure at the inlet and the outlet of the blood vessels. Finite
element method (FEM) is applied to solve the the Navier-Stokes (NS) equations.
The main components of the simulation can be divided into three parts: a)
dominant physics of the laminar fluid model, b) settings of the boundary
conditions, and c) development of a solver in space-time.

####  2.3.1 Laminar Physics

A three-dimensional NS equations with the assumption of stationary
incompressible fluid are applied in our study as the following:

$\rho(\boldsymbol{u}\cdot\nabla)\boldsymbol{u}=\nabla\left(-p\mathbb{I}+\mu[\nabla\boldsymbol{u}+(\nabla\boldsymbol{u})^{T}]\right)+\boldsymbol{F}$
(1)

$\rho\nabla\cdot\boldsymbol{u}=0$ (2)

where $\boldsymbol{u}$ is the flow velocity, $p$ is the pressure, $\mu$ is the
dynamic viscosity and $\rho$ the density of the blood, which is assumed to be
Newtonian fluid with $\mu$ = 0.005 Pa·s and $\rho$=1050 kg/m3[16].
$\boldsymbol{F}$ is the volume force vector and is set to be zero as we ignore
the gravity of blood. The Eq.1 shows the balance of momentum from Newton’s
Second Law, and the Eq.2 represents the constraint of continuity, which is
incompressible.

We use only laminar physics here rather than the turbulence and non-Newtonian
modelling, or Fluid-Structure Interaction (FSI) model, due to the small
Reynolds number (Re = 27) and our consideration of vascular sclerosis. On the
other hand, although the geometry of vessels counts a lot in the CFD model and
admittedly, it will be more accurate if we use FSI, patients suffered from
CSPH experience significantly reduced wall compliance, rendering our
assumption of rigid wall reasonable and computationally efficient.

####  2.3.2 Boundary Condition

The inlet is the entrance of the portal vein, with the velocity set as a
continuous function of time by interpolating cubic spline function into 30
uniformly-spaced time intervals in a cardiac cycle. Although the velocity has
already been averaged over time, we prolong the simulation time of interest by
duplicating the cardiac cycle four times (see in fig. 2 (b)), guaranteeing
that the flow is fully developed. A smooth transition from zero is added at
the beginning 100 ms to accelerate the convergence, as enlarged in fig. 2 (a).

The outlet is naturally the exit of the hepatic vein, and a reference pressure
is set here. As vHVPG is calculated as the pressure difference before and
after the balloon is inflated, it is irrelevant to the absolute value of the
reference pressure. For the sake of similarity to clinical measurement and
following the physiological condition, the constant value of 7.5 mmHg is set
as the outlet pressure condition.

Figure 2: Velocity settings in an example of a patient with liver cirrhosis.
(a) The illustration of how to smooth the abrupt transition at the beginning
100 ms of a flow velocity period. (b) The exact interpolated function of
velocity consists of four duplicate periods. (c) Corresponding vHVPG as a
function of time. The grey denotes the section enclosed by the function and
the x-axis, which is the integration of vHVPG during the calculation time
considering the plus and minus sign. Then the integration by the time is the
average vHVPG.

####  2.3.3 Solver

We use the Laminar Fluid module of COMSOL Multiphysics for finite element
simulation. We use a fully coupled direct solver for CFD models. To avoid run-
out-of-memory error, we cut the iteration to be ten thousand times. Two
dependent variables are coupled in the equation, i.e., the pressure and
velocity field. Moreover, algebraic multigrid (AMG) solvers are applied to
obtain these two values at each node[5].

###  2.4 Experiments

This retrospective study includes a total of ten healthy volunteers and five
patients with liver cirrhosis. The study was approved by the ethics committee
of the local hospital. Informed consent was obtained from all the subjects
prior to the examinations. The cirrhosis was diagnosis by means of liver
biopsy.

####  2.4.1 Building Geometry

Portal veins and hepatic veins are extracted from PC-MRI separately and
therefore the connection of capillary is needed. Since CFD largely depends on
the geometry, the strategy to rebuild it is fundamental and original here in
our work. As can be viewed in Fig.3, in principle, the connection of portal
vein and hepatic vein should get close to the reality, where both one-to-one
and one-to-many connections exist. Therefore, the decisions to carry out which
kind of the link and how to bridge the neighbour ports are important, and we
demonstrate our procedure as below.

Figure 3: Procedure to connect vessels by interpolation. The top line is the
case when the angle between two vessels $\theta>90\degree$, and the bottom
line is $\theta<90\degree$. $\overrightarrow{P_{0}P_{1}}$ defines the
directional vector $\boldsymbol{v_{1}}$ of the vessel on the left and
$\overrightarrow{P_{2}P_{3}}$ works for $\boldsymbol{v_{2}}$ on the right.
$P_{1}^{\prime}$ (or $P_{3}^{\prime}$) is the translation of $P_{1}$ (or
$P_{3}$) along $\boldsymbol{v_{1}}$ (or $\boldsymbol{v_{2}}$) at the distance
of $|P_{1}P_{3}|/2$. The midpoint $P_{c}$ is in the middle of $P_{1}^{\prime}$
and $P_{3}^{\prime}$. In (b) and (e), interpolation is demonstrated as green
dots along $P_{0},P_{1},P_{1}^{\prime},P_{c},P_{3}^{\prime},P_{3}$. As in (c),
large $\theta$ results in a smooth bridge and no more connections are allowed
for them, while in (f) the curve will provide a new point $P$, served as a new
port to be connected with a third vessel.

First of all, find two ends of the axis in two vessels, say
${P_{0},P_{1},P_{2},P_{3}}$, from which the direction vector
$\boldsymbol{v_{1}}$ and $\boldsymbol{v_{2}}$ will be obtained. Then defining
the center point $P_{c}$ by the midpoint of $P_{1}^{\prime}P_{3}^{\prime}$,
where the spot $P_{1}^{\prime}$ is translated from $P_{1}$ along the vector
$\boldsymbol{v_{1}}$ by the distance of $|P_{1}P_{3}|/2$ and the same for
$P_{3}^{\prime}$ but from $P_{3}$ along $\boldsymbol{v_{2}}$. Next,
interpolate enough points between the point set of
$P_{0},P_{1},P_{1}^{\prime},P_{c},P_{3}^{\prime},P_{3}$ and extrude the circle
surface with the diameter of initial diameter on both sides along the
interpolated curve. Finally, with some polishing, the connection is built up.
However, if the angle between two vectors $\boldsymbol{v_{1}}$ and
$\boldsymbol{v_{2}}$ is acute, the bridge will end up with a rather sharp
curve, generating a new point $P$ to be connected with other vessels.

Figure 4: Average pressure distribution of typical CFD model. The first two
pictures show the pseudo-color images of flow pressure for a healthy volunteer
and the other two are for a patient with liver cirrhosis. Arrows show
measurement location. As for virtual free hepatic venous pressure (vFHVP),
(a)=9.8 mmHg, (c)=9.6 mmHg; for virtual wedged hepatic venous pressure
(vWHVP), (b)=14.5 mmHg, (d)=19.8 mmHg. vHVPG = vWHVP - vFHVP, and it is 4.7
mmHg for the healthy person, but 10.2 mmHg for the patient with liver disease.
As CSPH is classified as HVPG $\geq$ 10 mmHg, the findings of these two cases
fulfill this standard, therefore verifying our method of CFD model based on
PC-MRI.

After the connection, a break point is created for the modeling of WHVP,
representing a strict block of the blood flow along the vessel. We ring cut
the right-most secondary branch at the point of about 1 cm away from the fork.
Afterwards, both sides are healed by the smooth smallest surface and then
follows the remesh process. Remeshing includes uniforming the surface mesh and
creating the volume mesh from the surface, providing an exact calculation
domain for CFD with FEM. See in Fig. 4, pressure is evaluated near the break
point from the side closer to the portal vein. The virtual FHVP (vFHVP) is
from the model of complete blood vessels, while the virtual WHVP (vWHVP) is
achieved in the cut-off model. Therefore, similar to the clinical measurement
in practice, the virtual hepatic venous pressure gradient (vHVPG) is the
difference between the above two pressures: vHVPG = vWHVP - vFHVP.

##  3 Results and Discussion

In this section we present the simulation result and prove the accuracy by
analyse the data. In fig. 4, typical examples of average pressure distribution
in a healthy volunteer and a patient with liver cirrhosis can be compared with
each other. Three main properties verify the CFD model primarily: 1) Flow
pressure decreases gradually from the inlet at the portal vein to the outlet
at the hepatic vein. 2) Averagely, vWHVP is larger than vFHVP in both cases.
3) Obviously, the flow pressure of the patient is larger than that of the
healthy volunteer. The simulation result of vHVPG is concluded in Table 1.

HVPG larger than 10 mmHg is defined as CSPH for the statistical analysis.
Although the number of cases is limited, it still shows a significant increase
of HVPG for the patients and 80 percents of them are found to be CSPH (HVPG
$\leq$ 10 mmHg) according to our studies. Based on this confusing matrix, the
kappa test has shown substantial consistency (k=0.7) between the diagnosis
based on histopathology and our virtual HVPG measurements. Besides, the normal
HVPG is suggested to be in the range of 3-5 mmHg for healthy adults[3], which
is also fulfilled in our studies.

Table 1: Summary of the vHVPG for Different Groups Group | Average vHVPG (mmHg) | No. of true | No. of false
---|---|---|---
Healthy volunteers | $5.9\pm 2.3$ | 9 ($\leq 10$mmHg) | 1 ($\geq 10$mmHg)
Patients with liver cirrhosis | $21.5\pm 8.8$ | 4 ($\geq 10$mmHg) | 1 ($\leq 10$mmHg)

As for the next step, the interaction between structures and the fluid will be
considered, resulting in a more precise model than the current one of simple
single-phase laminar physics. Also, in order to virtually model the the
invasive measurement as close as possible, the geometry ought to be improved
as well. Although the spatial resolution of CT is usually higher than MRI for
liver imaging, the high-resolution 3D GRE sequence used in this study achieves
resolution of 0.87 × 0.87 × 2.0 mm3. We believe that the resolution is
comparable to conventional CT and adequate for CFD modelling. Interestingly,
although the fine-scale capillaries are not modelled, the vHVPG result is
consistent with histopathology diagnosis. It has been reported that the side
branches do affect some specific wall variables of CFD simulation results,
e.g. wall shear stress, although it is in the cardiovascular system [8, 11,
15]. However, it was suggested that including side blood vessels up to 1 mm in
diameter is appropriate [15]. While the vessel diameter that can be observed
in our study is no less than 0.87 mm corresponding to the spatial resolution,
we believe such simplification has a limited influence on the calculation of
pressure gradient. However, we still believe that the way forward lies in some
modified CFD model including nearly all the vessels or at least their effects,
constructed automatically with as little manual interference as possible.

##  4 Conclusion

In this article, we introduced a one-stop solution to acquire an accurate
estimate of HVPG non-invasively, assisted by PC-MRI and CFD. For the first
time, CFD model is established entirely based on MRI for the geometry and
boundary conditions of velocity in the liver system, as we take advantage of
the latest technology of PC-MRI. Since the limit of resolution forbids a
complete map of entire blood vessels in the liver, trunks of hepatic vein and
portal vein are preserved and connected automatically by our original
programming method, which demonstrates complex networking very similar to the
real vessels. FEM solver is applied in the CFD model accelerating the
calculation speed. By now, it only takes less than twenty minutes to compute
HVPG from the raw image, while it can even be faster in our further study.
Testing in a data set of ten healthy people and five patients with liver
cirrhosis, we verify the accuracy of our method in diagnosing people with
CSPH. Our study facilitates the process of non-invasive and patient-specific
treatment.

####  4.0.1 Acknowledgements

This study was funded by the National Natural Science Foundation of China (No.
81971583, No. 62001120), National Key R & D Program of China (No.
2018YFC1312900, No. 2019YFA0709502), Shanghai Natural Science Foundation (No.
20ZR1406400), Shanghai Municipal Science and Technology Major Project
(No.2017SHZDZX01, No.2018SHZDZX01), Shanghai Municipal Science and Technology
(No.17411953600), Shanghai Sailing Program (No. 20YF1402400), ZJLab and Key
Laboratory of Computational Neuroscience and BrainInspired Intelligence (Fudan
University), Ministry of Education, China.

## References

  * [1] Castéra, L., Bail, B.L., Roudot-Thoraval, F., Bernard, P.H., Foucher, J., Merrouche, W., Couzigou, P., Lédinghen, V.D.: Early detection in routine clinical practice of cirrhosis and oesophageal varices in chronic hepatitis c: comparison of transient elastography (fibroscan) with standard laboratory tests and non-invasive scores. Journal of Hepatology 50(1), 59–68 (2009)
  * [2] Cito, S., Pallarés, J., Vernet, A.: Sensitivity analysis of the boundary conditions in simulations of the flow in an aortic coarctation under rest and stress conditions. In: Camara, O., Mansi, T., Pop, M., Rhode, K., Sermesant, M., Young, A. (eds.) Statistical Atlases and Computational Models of the Heart. Imaging and Modelling Challenges. pp. 74–82. Springer Berlin Heidelberg, Berlin, Heidelberg (2014)
  * [3] De, F.R., Faculty, B.V.: Expanding consensus in portal hypertension: Report of the baveno vi consensus workshop: Stratifying risk and individualizing care for portal hypertension. Journal of Hepatology 63(3), 743–752 (2015)
  * [4] Dong, J., Qi, X.: Liver imaging in precision medicine. Ebiomedicine p. S2352396418301816 (2018)
  * [5] Falgout, R.D.: An introduction to algebraic multigrid computing. Computing in Science and Engineering 8(6), 24–33 (2006)
  * [6] Fang, C., An, J., Bruno, A., Cai, X., Qi, X.: Consensus recommendations of three-dimensional visualization for diagnosis and management of liver diseases. Hepatology International 14(1) (2020)
  * [7] Franchis, R.D.: Portal Hypertension VI: Proceedings of the Sixth Baveno Consensus Workshop: Stratifying Risk and Individualizing Care. Springer (2016)
  * [8] Giannopoulos, A.A., Chatzizisis, Y.S., Maurovich-Horvat, P., Antoniadis, A.P., Hoffmann, U., Steigner, M.L., Rybicki, F.J., Mitsouras, D.: Quantifying the effect of side branches in endothelial shear stress estimates. Atherosclerosis 251, 213—218 (August 2016)
  * [9] Gorelick, P.B., Wong, K.S., Bae, H.J., Pandey, D.K.: Large artery intracranial occlusive disease: a large worldwide burden but a relatively neglected frontier. Stroke; a journal of cerebral circulation 39(8), 2396 (2008)
  * [10] Jenkins, R.H., Mahal, R., Maceneaney, P.M.: Noninvasive imaging of carotid artery disease: critically appraised topic. Canadian Association of Radiologists Journal 54(2), 121–123 (2003)
  * [11] Li, Y., Gutiérrez-Chico, J.L., Holm, N.R., Yang, W., Hebsgaard, L., Christiansen, E.H., Mæng, M., Lassen, J.F., Yan, F., Reiber, J.H., Tu, S.: Impact of side branch modeling on computation of endothelial shear stress in coronary artery disease: Coronary tree reconstruction by fusion of 3d angiography and oct. Journal of the American College of Cardiology 66(2), 125–135 (2015)
  * [12] Liu, F., Ning, Z., Liu, Y., Liu, D., Tian, J., Luo, H., An, W., Huang, Y., Zou, J., Liu, C.a.: Development and validation of a radiomics signature for clinically significant portal hypertension in cirrhosis (chess1701): a prospective multicenter study. Ebiomedicine (2018)
  * [13] Schaller, J., Goubergrits, L., Yevtushenko, P., Kertzscher, U., Riesenkampff, E., Kuehne, T.: Hemodynamic in aortic coarctation using mri-based inflow condition. In: Camara, O., Mansi, T., Pop, M., Rhode, K., Sermesant, M., Young, A. (eds.) Statistical Atlases and Computational Models of the Heart. Imaging and Modelling Challenges. pp. 65–73. Springer Berlin Heidelberg, Berlin, Heidelberg (2014)
  * [14] Tana, M.M., Muir, A.J.: Diagnosing liver fibrosis and cirrhosis: Serum, imaging, or tissue? Clinical Gastroenterology and Hepatology (2018)
  * [15] Vardhan, M., Gounley, J., Chen, S.J., Kahn, A.M., Leopold, J.A., Randles, A.: The importance of side branches in modeling 3d hemodynamics from angiograms for patients with coronary artery disease. Scientific reports 9(1), 8854 (2019)
  * [16] Wadhawan, M., Dubey, S., Sharma, B.C., Sarin, S.K., Sarin, S.K.: Hepatic venous pressure gradient in cirrhosis: Correlation with the size of varices, bleeding, ascites, and child’s status. Digestive Diseases and Sciences 51(12), 2264–9 (2006)
  * [17] Zhang, X., Zhao, X., Tan, T.: Robust dialog state tracker with contextual-feature augmentation. Applied Intelligence 51(2-3), 2377–2392 (2021)
  * [18] Zhao, X., Liu, Y., Ding, J., Zhang, M., Fu, W., Bai, F., Ren, X., Qiao, A.: Multiscale study on hemodynamics in patient-specific thoracic aortic coarctation. In: Camara, O., Mansi, T., Pop, M., Rhode, K., Sermesant, M., Young, A. (eds.) Statistical Atlases and Computational Models of the Heart. Imaging and Modelling Challenges. pp. 57–64. Springer Berlin Heidelberg, Berlin, Heidelberg (2014)