Sofiane Aboura
email: sofiane.aboura@univ-paris13.fr

Emmanuel Lépinette
email: emmanuel.lepinette@ceremade.dauphine.fr

Evaluation of the Fair Credit Risk Premium in Commercial Lending

Keywords: Corporate Loans, Credit Risk, Investment Decisions, Portfolio Choice, Contingent Pricing 2000 MSC: 60G44, G11-G13, G32-G33

We consider the problem of characterizing and computing the fair credit risk premium that a firm should pay when borrowing money from a bank. Using a risk-neutral approach, we show that there is a unique credit risk premium for a commercial loan depending on the firm's strategy such that the expected discounted value of the bank's payoff coincides with the loan's par value. We then propose a numerical procedure to estimate the premium.

Introduction

This paper discusses the pricing of corporate loans in the case where the risk premium is explicitly related to the firm's investment strategy.

A corporate loan is seen as a bilateral financial agreement issued, in the form of periodic payments, by a bank for its corporate non-financial firms. The main difference between loans and bonds is that banks are senior to bondholders in bankruptcy [START_REF] Schwert | Is Borrowing from Banks More Expensive than Borrowing from the Market? Working Paper[END_REF]. Loans are typically unsecured, which places the bank on equal level with unsecured bondholders, although it is unusual for a firm to default on an unsecured loan 1 . Of course, a firm can also borrow money by issuing unsecured corporate bonds in the loan market. To that end, [START_REF] Santos | Bank Corporate Loan Pricing Following the Subprime Crisis[END_REF] finds that after the sub-prime crisis, banks increased the interest rates on their loans to bank-dependent borrowers by more than they did on their loans to borrowers that have access to the bond market. This result is also confirmed by [START_REF] Schwert | Is Borrowing from Banks More Expensive than Borrowing from the Market? Working Paper[END_REF] who shows that firms are willing to pay a higher cost to borrow from a bank due to the bank specialness. In clear, it seems there exists a preference for banks despite a relatively higher cost. The overall cost for a loan can be seen as the basis point spread over the money market rate (the main reference rate for corporate loans is Libor) plus the annual fees over the life of the loan.

The pricing of loan contracts has attracted considerable research interest 2 . The pricing structure can depend on the type of loans (Credit Lines versus Term Loans) as explained by [START_REF] Berg | Mind the Gap: The Difference between U.S. and European Loan Rates[END_REF] 3 as some loan contracts are set with repayment on maturity (zero-coupon plan), while others are designed with amortization plans (periodical payments plan).

An important problem with loan pricing comes from the assessment of the default probability. Credit risk is generally calculated implicitly from the price of loans, as it includes the investors expectation about future possible losses. But in practice, this information is provided by credit-rating agencies that compute long-term historical default rates for various horizons of time. These corporate credit ratings reduce the gap in terms of information asymmetry between lenders and borrowers to find an equilibrium price [START_REF] Blöchlinger | Credit Rating and Pricing: Poles Apart[END_REF]. The corresponding credit scoring is an important instrument in setting an appropriate default premium when determining the rate of interest charged to the borrower company with a past credit history. Hence, the loan spreads reflect the ex-ante credit risk as they are strongly correlated with probabilities of borrower default and capture the loan's credit risk smoothed through the cycle [START_REF] Lee | Risk Taking and Interest Rates: Evidence from Decades in the Global Syndicated Loan Market[END_REF]. However, the credit rating does not incorporate the forward looking market consensus of the default risk and overall, it seems that past corporate bond spreads appear to have poor predictive power over default rates [START_REF] Giesecke | Corporate bond default risk: A 150-year perspective[END_REF]. A related prominent problem with loan contracts is the determination of the fair risk premium that reflects the level of risk of the borrower's project. Each risky loan carries a risk premium for which the current price is the expected value of the pay-off discounted at the risk-free interest rate augmented by a risk premium. Loan risk premiums correspond to the loan rates less the risk-free rates and lower risk premiums arise when they carry lower losses due to default. In case of high risk profile for the borrowing firm, the bank can limit its exposure by requiring the riskier borrowers to select only short-term loans or imposing a collateral asset for hedging purpose. Interestingly, Allen, Scott and Vasso (2016) explain that in fact lenders offer a menu of contract terms such that applicants with higher-quality projects choose secured debt with lower risk premiums, while applicants with lower-quality projects selfselect into unsecured debt with higher risk premiums. The valuation problem becomes relying on the firm's investment strategy given the relation between the quality of a project and the associated risk premium.

To that end, the present paper considers credit risk from the perspective of a bank willing to lend money according to the strategy of a non-financial firm. The very question is how large the minimum credit interest rate the bank re-quires should be for it to compensate for the risk not to be reimbursed given the firm's strategy? This paper considers the realistic case where the bank has to fix the fair risk premium depending on the firm's strategy expressed in terms of quality of investment and level of consumption.

The article is organized as follows. Section 2 displays the literature review. Section 3 outlines the model. In Section 4, the main results are presented with their mathematical proofs. Section 5 concludes.

Review of literature

In financial theory, credit risk is defined as the case where the asset value ends up below the level of debts driving the firm to be insolvent. The fundamental contribution of the Merton (1974)'s structural model is that firm' risky debt can be viewed as a pseudo-bond equal to a portfolio composed of a risk-free debt minus a short put option on the value of the firm' assets. Structural models jointly capture the capital structure along with the evolution of firm asset value [START_REF] Lee | Risk Taking and Interest Rates: Evidence from Decades in the Global Syndicated Loan Market[END_REF]. The structural model gives rational restrictions on the pricing of corporate debt and equity claims on a levered firm through classical option pricing theory (Black and Scholes, 1973;Merton, 1973) and classical corporate finance theory [START_REF] Modigliani | The cost of capital, corporation finance and the theory of investment[END_REF]. This derivatives-based approach gives rational for pricing loans under a risk-neutral framework.

Important extensions of the Merton (1974)'s structural model has been done so far (Geske, 1977;Ingersoll, 1977) and particularly for credit risk because Merton model credit spreads are too low regarding what market prices. The literature refers to this phenomenon as the "credit spread puzzle" to explain why actual credit spreads are well above the theoretical spreads predicted by the debt valuation models begining from the seminal Merton (1974) model.

There is an extensive theoretical literature examining this credit spread puzzle through the pricing of bonds in the presence of credit risk. Among the main pieces of work, study of the effects of safety covenants giving the bondholders the right to bankrupt and restrictions on the financing of interest and dividend payments (Black and Cox, 1976); pricing of contingent claims on the term structure of interest rates (Ho and Lee, 1986; Heath, Jarrow, Morton, 1992); derivation of solution for risky debts for finite maturity and stochastic risk-free interest rate (Longstaff and Schwartz, 1992); connection of debt value to firm risk, taxes, bankruptcy cost and bond covenants (Leland, 1994); valuation of risky corporate debt that incorporates both default and interest rate risk (Longstaff and Schwartz, 1995); pricing of derivatives securities involving credit risk (Jarrow and Turbull, 1995); optimal default [START_REF] Leland | Optimal Capital Structure, Endogenous Bankruptcy, and the Term Structure of Credit Spreads[END_REF]; pricing of interest-rate swap yields with a standard term structure model incorporating default risk (Duffie and Singleton, 1997); agency costs [START_REF] Leland | Agency Costs, Risk Management, and Capital Structure[END_REF]; asymmetric information and uncertainty [START_REF] Duffie | Term Structures of Credit Spreads with Incomplete Accounting Information[END_REF]; relation between risk premium in corporate bond prices and firms' investment and financing policies (Kuehn and Schmid, 2014); complex debt structures with multiple bonds and various covenants (Liu, Dai and Wang, 2016); idiosyncratic asset uncertainty measure as the residual volatility from a market model on equity [START_REF] Culp | Option-Based Credit Spreads[END_REF].

However, Collin-Dufresne, Goldstein, and Martin (2001) find that credit spreads are driven by factors difficult to explain by standard credit models, which has been confirmed by Giesecke, Longstaff, Schaefer, and Strebulaev (2011) who study corporate bond default rates using an extensive data set (1866-2008). In fact, several assumptions behind these structural models tend to deviate the results from real world observations [START_REF] Liu | Evaluating corporate bonds and analyzing claim holders' decisions with complex debt structure[END_REF]. For instance, the structural models tend to simplify the capital structure to build mathematical models that are tractable but fail to capture the complexity of debt structure [START_REF] Eom | Structural models of corporate bond pricing: An empirical analysis[END_REF]). In addition, structural models have difficulty in simultaneously explaining both the default probabilities and true credit spreads particularly for investment-grade bonds [START_REF] Huang | How much of the corporate-treasury yield spread is due to credit risk?[END_REF]. This paper adds to the broad literature on structural models of firm by deriving the fair value of the credit risk premium as being explicitly related to the firm's strategy. The novelty is to define the firm's strategy as a combination of the investment strategy (portfolio of risky assets) along with the financial policy (deleveraging ratio) and the dividend policy (pay-out ratio). Whatsoever, this fair pricing relies on the ability of the bank to share the same level of information with the borrowing non-financial firm.

The model

Let T > 0 be an horizon date. We consider a financial market model defined on a complete stochastic basis (Ω, (F t ) t∈[0,T ] , P) satisfying the usual assumptions. Let (r 0 t ) t∈[0,T ] be a stochastic process adapted to (F t ) t∈[0,T ] interpreted as the risk-less interest rate of the market, i.e. such that it is possible for any firm to borrow money on the financial market at the interest rate r 0 if and only if the firm is in position to reimburse its debt with probability 1. We denote by S 0 the risk-less asset S 0 defined by S 0 0 = 1 and the dynamics dS 0 t = r 0 t S 0 t dt, i.e.:

S 0 t = exp t 0 r 0 u du , t ∈ [0, T ]. (3.1)
Let us consider a firm characterised (as in the Merton model) by its debt (D t ) t∈[0,T ] and the asset (A t ) t∈[0,T ] so that the equity is given by (E t ) t∈[0,T ] such that E = (A-D) + . We suppose that the stochastic process D is adapted to (F t ) t∈[0,T ] and it is defined by its initial value D 0 , i.e. the amount of capital the firm needs to borrow at time 0, and the stochastic differential equation

dD t = r t D t dt -K t dt, t ∈ [0, T ]. (3.2)
where r t is the instantaneous interest rate of the credit at time t and K t is the instantaneous amount of money that the firm reimburses at time t per unit of time. In the following, we set K t := k t D t so that k t is the ratio of the debt D t with is reimbursed at time t.

The asset A of the firm satisfies by definition A t = A (θ 0 ,θ) t := θ 0 t S 0 t + θ t S t where θ 0 is the quantity of asset S 0 held by the firm. Similarly, θ t is the vector of all quantities of risky assets S = (S 1 , • • • , S d ) invested by the firm at time t. Notice that the quantity θ t S t is the usual scalar product between the vectors θ t and S t . We suppose that S is a semi-martingale and the processes θ 0 and θ are predictable. Moreover, θ is supposed to be integrable with respect to S. In the following, we suppose the self-financing condition

dA t = θ 0 t dS 0 t + θ t dS t -c t dt, t ∈ [0, T ], (3.3)
where c t ≥ K t is a stochastic process adapted to (F t ) t∈[0,T ] . Notice that Ãt = Ẽt + Dt is supposed to be positive. We interpret the difference c t -K t ≥ 0 as the amount of dividends delivered to equity holders at time t per unit of time. The financial interpretation of (3.3) is clear: there is no adding or withdrawn of wealth but the amount of debt which is reimbursed plus the amount of dividends paid on the period [t, t + dt] are withdrawn from the firm's asset.

In the following, for any (F t ) t∈[0,T ] -adapted stochastic process X, we use the notation Xt = X t (S 0 t ) -1 for the discounted value of X. As (S 0 ) -1 is of finite variations and continuous, precisely d(S 0 t ) -1 = -r 0 t (S 0 t ) -1 dt, we deduce by the integration by parts formula that

d Xt = -r 0 t Xt dt + (S 0 t ) -1 dX t , t ∈ [0, T ]. (3.4)
For our model, we get that

d Ãt = θ t d St -ct dt, t ∈ [0, T ], (3.5
) d Dt = (p t -k t ) Dt dt, t ∈ [0, T ], (3.6)
where p t := r t -r 0 t is the credit risk premium. When p = 0, the debt process D is denoted by D 0 . In the following, we only consider admissible strategies θ such that Ãt ≥ -κ θ for all t ∈ [t, T ] where κ θ ≥ 0 is a constant. Observe that, in practice, κ θ is generally equal to 0. Moreover, contrarily to the Merton model where the asset A is defined ex ante, the dynamics above shows that A may depend on the risk premium through the process c ≥ k D. This is very natural since we may think that the debt interest should adversely impact the firm's asset. Moreover, the firm is characterised by the financial position (-D, Ã) whose liquidation value at time t is Lt := Ãt -Dt so that Ẽ = L+ . At time 0, we suppose that Ã0 -D0 ≥ 0. Observe that the dynamics of L is given by

d Lt = θ t d St -dt dt -p t Dt dt, t ∈ [0, T ], (3.7)
where dt = ct -k t Dt is the amount of dividends given to the share holders on the period [t, t + dt]. The dynamics above shows that the liquidation value of the firm's financial position is naturally controlled by the investment strategy θ but it is adversely impacted by the dividends d ≥ 0 paid to the share holders and by the credit risk premium p as well. In particular, apart the risk provided by the risky asset S, there is a risk generated by the strategy θ such that an increase of the credit risk premium may decrease the liquidation value L, which may leads to a bankruptcy when L = 0.

Observe that, by In the case where the inequality above holds, then it also holds for p = 0. Moreover, as D is an increasing function of p, the inequality is violated as soon as p is large enough. Therefore, it is possible for the debt holders to deliberately make the firm insolvent when increasing the credit risk premium. It is natural to use the convention that the credit risk premium is fixed to 0 as soon as the left hand side of the inequality above is non negative for p = 0. This convention is actually justified under the condition (MP) below, see Lemmas 4.4 and 4.9.

Fair credit risk premium

We define the fair credit risk premium as the smallest process p = p such that the inequality (3.8) holds, i.e. ÃT ≥ DT where D satisfies the dynamics (3.6) with p = p. As Ã0 = Ẽ0 + D0 , the conclusion follows. Reciprocally, if (4.9) holds, then ÃT ≥ DT when p = 0. 2

Application to a complete market model

Suppose that there exists Q ∼ P such that S is a martingale under Q. We suppose that the market is complete, i.e. any integrable contingent claim may be replicated by a portfolio process which is a martingale under Q. Therefore, it is well known that Q is the unique risk-neutral probability measure for S. The dividend plans d are supposed to be admissible in the sense that the cumulated dividend plan T 0 dt dt is supposed to be integrable under Q. Inspired by Lemma 4.1, since the market is complete, consider E * 0 and the admissible investment strategy θ * such that

E * 0 + t 0 θ * u d Su = E Q T 0 du du|F t , t ∈ [0, T ], a.s. (4.10)
Observe that, if E 0 = E * 0 , then θ * is a strategy which allows a zero credit risk premium p = 0 by Lemma 4.1. The following proposition shows that θ * is the unique strategy that yields the zero credit premium p = 0 when E 0 = E * 0 . Proposition 4.2. Suppose that E 0 = E * 0 . Then, ÃT ≥ DT if and only if θ = θ * d[ S, S] a.s. In particular, the investment strategies θ and θ * generate the same asset Ã.

Proof. Suppose that ÃT ≥ DT and E 0 = E * 0 . Using (3.8) and (4.10), we deduce that

T 0 (θ u -θ * u )d Su ≥ 0, a.s.
Since θ is admissible, and c ≥ 0, we deduce the existence of a constant κ ≥ 0 such that

M t := t 0 (θ u -θ * u )d Su ≥ -κ - t 0 θ * u d Su , t ∈ [0, T ].
By (4.10), the process in the r.h.s. of the inequality above is a martingale. Therefore, using the Fatou's lemma, we deduce that the local martingale M is a supermartingale. Since M T ≥ 0, we deduce that M t ≥ 0 for every t ∈ [0, T ]. On the other hand, M 0 = 0 and E Q (M t |F 0 ) ≤ M 0 by the super martingale property. Therefore, M t = 0 a.s. for all t ≤ T and, by [10, 4.54],

0 = E Q M 2 T = E Q T 0 (θ t -θ * t ) 2 d[ S, S] t .
We deduce that θ = θ * d[ S, S] a.s. Reciprocally, if this property holds, then M = 0 hence (3.8) holds with p = 0. 2

Following the same arguments, we obtain the following general result:

Proposition 4.3. ÃT ≥ DT if and only if E 0 ≥ E * 0 and E 0 + t 0 θ u d Su ≥ E * 0 + t 0 θ * u d Su , ∀t ∈ [0, T ].
Proof. As θ is supposed to be an admissible strategy, we deduce that the local martingale t → t 0 θ u d Su is a supermartingale. Therefore, with the inequality ÃT ≥ DT , i.e. (3.8), we obtain that:

E 0 + t 0 θ u d Su ≥ E 0 + E Q T 0 θ u d Su |F t , ≥ E Q T 0 du du|F t , ≥ E * 0 + t 0 θ * u d Su .
With t = 0, we deduce that E 0 ≥ E * 0 . This proposition above means that the debt holder should require a non null credit risk premium as soon as the investment strategy is not efficient enough or when the initial equity capital is too small.

Valuation under the risk-neutral probability measure

In the following, we suppose that A ≥ 0. Let Q ∼ P be the risk-neutral probability measure for S in the market model supposed to be complete. In the case of a possible default, e.g. when p = 0, the discounted payoff delivered to the debt holders is

hD T := T 0 k t Dt dt + ÃT ∧ DT , (4.11)
i.e. the sum of the cumulated amount of money reimbursed before T by the firm and the residual debt DT when it is smaller than ÃT . Otherwise, in the case of bankruptcy, the debt holders receive the asset value ÃT . The payoff delivered to the equity holders is

hE T := T 0 du du + ( ÃT -DT ) + , (4.12)
i.e. the sum of the cumulated dividends and the asset value minus the debt to be reimbursed. Observe that:

hD T + hE T = E 0 + D 0 + T 0 θu d Su , (4.13)
i.e. the terminal gain generated by the investment strategy θ when starting from the initial endowment E 0 +D 0 (i.e. the r.h.s. of (4.13)) is shared between the debt holders and the equity holders with a priority to the debt holders, see (4.11). Since the strategy θ is supposed to be admissible, the

Q-local martingale t → t 0 θ u d Su is a Q-supermartingale hence E Q t 0 θ u d Su ≤ 0 for every t ≤ T .
It is natural to require a stronger condition, i.e. E Q T 0 θ u d Su = 0 since the investors are willing to optimize their profits. Therefore, under this condition, N : Proof. Using (4.13), we deduce that E Q ( hD Recall that x -= -min(x, 0). Notice that hD T (0) ≤ D 0 so that the debt holders should accept to lend money at the rate p = 0 if and only if hD Remark 4.6. Let us suppose that the inequality hD T (0) = D0 is not satisfied a.s. Therefore, P ( hD T (0) < D 0 ) > 0 and hD T (0) ≤ D 0 a.s. implies that E Q ( hD T (0)) < D 0 . Therefore, it is necessary to increase the credit risk premium p for Condition (MP) to be satisfied.

t → t 0 θ u d Su is a martingale under Q. Indeed, if E Q (N T |F t ) ≤ N t and if the inequality is strict on a non null set, we get that E Q (N T ) < E Q (N t ) ≤ 0. So, the condition E Q (N T ) = 0 yields a contradiction.
T ) + E Q ( hE T ) = D 0 + E 0 . This equality implies that E Q ( hD T ) = D 0 and E Q ( hE T ) = E 0 as soon as E Q ( hD T ) = D 0 . 2 Recall that Ãt = A 0 + t 0 θ u d Su - t 0 du du - t 0 k u Du du, Dt = D 0 + t 0 (p u -k u ) Du du. It follows that ÃT ≥ DT ⇔ A 0 + T 0 θ u d Su - T 0 p u Du du ≥ D 0 + T 0 du du, ⇔ E 0 + T 0 θ u d Su -
T (0) = D 0 , i.e. if A 0 + T 0 θ u d Su -

Constant credit risk premium

In this section, we only consider constant credit risk premiums p u = p 0 = p for every t ∈ [0, T ]. The Merton model is a particular case where k = d = 0 and the discounted face value KT of the debt is a given deterministic constant. We then deduce the unique credit risk premium p by solving the equation KT = D 0 e pT where D 0 is the market value of hD T = ÃT ∧ KT . In the following, we show the uniqueness of such a constant credit risk premium for our model in the case where E Q T 0 dr dt < E 0 . Notice that, if this inequality does not hold, there is no economical reason for the debt holders to lend money as the firm is not solvent (see Lemma below) even if the condition (MP) actually holds whatever p. In that case, we take the convention that the credit risk prime is p * = 0 as hD T (p) does not depend on p.

Proposition 4.7. Suppose that

E Q T 0 dr dt ≥ E 0 . Then, E Q ( hD T (p)) ≥ D 0 if and only if E Q T 0 dr dt = E 0 , and LT = ÃT -DT ≤ 0 a.s. i.e. E 0 + T 0 θ u d Su - T 0 du du ≤ T 0 p u Du du, a.s.

Under these conditions, hD

T (p) does not depend on p and E Q ( hD

T (p)) = D 0 . Proof. By (4.14), E Q ( hD T (p)) ≥ D 0 if and only if E Q (( LT ) -) ≤ E Q T 0 p u Du du or equivalently if E Q T 0 p u Du du -γ + ≤ E Q T 0 p u Du du ,
where Note that the function δ(x) = x -(x -γ) + is non decreasing and for x ≥ 0, |δ(x)| ≤ |γ|. Supposing that γ is integrable, we deduce by the dominated convergence theorem that φ(∞) = E Q (γ). In particular, since

γ = L 0 + T 0 θ u d Su - T 0 du du = LT + T 0 p u Du du. Since x + ≥ x, we deduce that E Q T 0 p u Du du -γ + ≥ E Q T 0 p u Du du -E Q (γ). Since L 0 = E 0 and E Q T 0 dr dt ≥ E 0 , then E Q (γ) ≤ 0. It follows that E Q ( hD T (p)) ≥ D 0 if and only if E Q (γ) = 0, i.e. E Q T 0 dr dt = E 0 , and
E Q T 0 p u Du du -γ + = E Q T 0 p u Du du . Therefore, E Q ( hD T (p)) ≥ D 0 if and only if E Q (γ) = 0 and
E Q (γ) = E 0 -E Q T 0 dt dt , we get that φ(∞) > 0. Moreover, φ(0) = -E Q γ -≤ 0. Therefore, as p → E Q ( hD T (p)
) is continuous, there exists p * ∈ R + such that φ(p * ) = 0. Suppose that there are two real numbers p 1 , p 2 ∈ R + such that φ(p 1 ) = φ(p 2 ) = 0 and p 1 < p 2 . Since δ is strictly increasing on (-∞, γ) and constant on [γ, ∞), we obtain that

δ T 0 p 1 Du (p 1 )du ≤ δ T 0 p 2 Du (p 2 )du
and, finally, the equality holds due to the assumption. Therefore, we necessarily have

T 0 p 2 Du (p 2 )du ≥ γ
since δ is strictly increasing on (-∞, γ). We deduce that φ(p 2 ) ≥ E Q (γ) where E Q (γ) > 0 by assumption. This yields a contradiction. 2

Note that, when φ(0) = 0, then γ = LT ≥ 0 a.s. hence ÃT ≥ DT a.s. so that, by the proposition below, p * = 0 is the only risk premium satisfying φ(p * ) = 0 under (MP). Note that g is a strictly increasing function on [0, ∞) such that g(0) = 0 and g(∞) = ∞. As g is also continuous, the inverse g -1 exists and is [0, ∞)valued.

We suppose that E Q T 0 dr dt < E 0 hence there is a unique credit risk prime p * ∈ R + by Proposition 4.8 such that φ(p * ) = 0. Note that p * coincides with the unique solution to the equation g(p) = E Q (γ -g(p)) -where γ is defined by

γ := E 0 + T 0 θ t d St - T 0 dt dt.
Here, we suppose that γ is square integrable. Let us consider the sequence p n defined by p 0 = 0 and

p n+1 = g -1 E Q (γ -g(p n )) -,
i.e. g(p n+1 ) = E Q (γ -g(p n )) -. As p 0 = 0, we have p 1 ≥ p 0 . Since g is increasing and x → x -is non decreasing, we may show by induction that p n+1 ≥ p n for all n ≥ 0. Similarly, as p 0 ≤ p * , we deduce that g(p 1 ) ≤ E Q (γ -g(p * )) -= g(p * ) hence p 1 ≤ p * and, by induction, p n ≤ p * for all n ≥ 0. We finally deduce that the increasing sequence (p n ) n converges to p * .

In practice, it seems to be difficult to implement the iteration scheme p n+1 = g -1 (E Q (γ -g(p n )) -) as we need first to estimate the expectation of (γ -g(p n )) -. Note that, since g is invertible, the problem is equivalent to solve the equation x = E Q (γ -x) -whose solution is x * = g(p * ) so that p * = g -1 (x * ). Clearly x * does not depend on D 0 while p * does as it is the solution to the equation x * = g(p * ).

In order to estimate x * , we propose a method inspired by the well known stochastic gradient descent. Precisely, we consider a sequence (γ n ) n≥1 of independent random variables distributed as γ and we introduce the following scheme: x 1 = 0 and

x k+1 = x k -α k H(γ k , x k ), k ≥ 0, where H(v, x) = x -(v -x) -,
(x) = EH(γ, x). We have C(0) = -E Q γ -and C(∞) = E Q γ >
0 by assumption. Recall also that C(x * ) = 0 by assumption. Since γ is integrable, we also deduce that C is continuous. Therefore, as C is non decreasing and the equation C(x) = 0 admits a unique solution, we deduce that, for every > 0,

inf |x-x * |> (x -x * )C(x) > 0. (4.16)
On the other hand, let us introduce D(x) = E Q H 2 (γ, x). We have

D(x) = E Q (H(x, γ) -H(x * , γ) + H(x * , γ)) 2 , ≤ 2E Q (H(x, γ) -H(x * , γ)) 2 + 2E Q H 2 (x * , γ), ≤ 8(x -x * ) 2 + 2E Q γ 2 .
The last inequality is deduced from the fact that

H(•, γ) is 2-Lipschitz and that |H(x, γ)| ≤ |γ| if x ≥ 0. Let us introduce the sequence d k = (x k -x * ) 2 , k ≥ 1 and the filtration F k = σ(γ 1 , • • • , γ k-1 ), k ≥ 2,
and F 1 is the trivial σ-algebra induced by the negligible sets. Notice that x k ∈ L 0 (R, F k ) for all k. The goal is to show that d k converges a.s. to zero. To do so, let us first compute

d k+1 -d k = (x k+1 -x * ) 2 -(x k -x * ) 2 = (x k -α k H(γ k , x k ) -x * ) 2 -(x k -x * ) 2 , d k+1 -d k = α 2 k H 2 (γ k , x k ) -2α k H(γ k , x k )(x k -x * ).
We deduce that

E(d k+1 -d k |F k ) = α 2 k E Q H 2 (γ k , x k )|F k -2α k (x k -x * )E Q (H(γ k , x k )|F k ) .
As γ k is independent of F k and x k is F k -measurable, we deduce that

E Q H 2 (γ k , x k )|F k = E Q H 2 (γ, x k ) , E Q (H(γ k , x k )|F k ) = E Q (H(γ, x k )) .
Therefore, by the first step,

E Q (d k+1 -d k |F k ) ≤ 8α 2 k (x k -x * ) 2 + 2α 2 k E Q γ 2 -2α k (x k -x * )E Q (H(γ, x k )) . (4.17)
Therefore, as α k (x k -x * )E Q (H(γ, x k )) ≥ 0,

E Q (d k+1 -d k (1 + 8α 2 k )|F k ) ≤ 2α 2 k E Q γ 2 . (4.18)
Let us introduce δ k+1 = µ k d k+1 , k ≥ 0, where µ k := Π k-1 j=1 (1+8α 2 k ) -1 , k ≥ 1 and µ 0 = 1. Multiplying (4.18) by µ k , we obtain that

E Q (δ k+1 -δ k |F k ) ≤ 2α 2 k µ k E Q γ 2 . (4.19)
From above, we deduce that the process M k := δ k -2E Q γ 2 k j=1 α 2 k µ k , k ≥ 1, is a super martingale. As δ k ≥ 0 for all k ≥ 1, M is bounded from above hence (M k ) k≥1 is a.s. convergent by the Doob's theorem. We deduce that (δ k ) k≥1 is also convergent. As the sequence (ln(µ k )) k≥1 is convergent, we deduce that µ k → µ > 0 hence the sequence d k = (x k -x * ) 2 is convergent.

Taking the expectation in (4.17) and summing up, as d k ≥ 0, we get that

E Q ∞ k=1 α k (x k -x * )E Q (H(γ, x k )) < ∞,
and finally

0 ≤ ∞ k=1 α k (x k -x * )E Q (H(γ, x k )) < ∞,
where we recall that α k (x k -x * )E Q (H(γ, x k )) ≥ 0 a.s. for all k ≥ 1.

Let us define r = lim inf k→∞ (x k -x * )E Q (H(γ, x k )). Suppose that r > 0. Then, (x k -x * )E Q (H(γ, x k )) > r/2 for k ≥ k 0 large enough so that

∞ k=k 0 α k (x k -x * )E Q (H(γ, x k )) > r/2 ∞ k=k 0 α k = +∞,
which yields a contradiction. We deduce that r = 0 and, in the case where lim k→∞ d k = 2 2 > 0, we get that (x k -x * ) 2 ≥ 2 for k large enough. By (4.16), we deduce a constant c such that (x k -x * )C(x k ) > c where we recall that C(x k ) = E Q (H(γ, x k )) hence a contradiction, i.e. lim k→∞ d k = 0. 2

Conclusion

We consider the problem of characterizing and computing the fair credit risk premium that a non-financial firm should pay when borrowing money from a bank. In particular, we study the case where the risk premium depends explicitly on the firm's strategy such that the expected discounted value of the bank's payoff coincides with the loan's par value. The strategy is explicitly defined in terms of investment and consumption. We use a riskneutral framework to show the existence of a unique credit risk premium for such a commercial loan. We then propose a numerical procedure to estimate the fair premium. A limit of this research is that it relies on the assumption of the absence of information asymmetry between the bank and the nonfinancial firm, which suggests a new direction for further research.

Lemma 4 . 1 .0 0 kt

 410 Let us consider a firm (D, A) as in Section 3. The fair credit risk premium is zero if and only if Ẽ0 + T 0 θ t d St ≥ T Dt dt, a.s. (4.9) Proof. If the credit risk premium is p = 0, then by (3.6), DT = D0 -T Dt dt. Moreover, by definition, p = 0 implies that ÃT ≥ DT . By (3.5)Dt dt.

Lemma 4 . 4 .

 44 Suppose that N : t → t 0 θ u d Su is a martingale under Q. Then, E Q ( hD T (p)) = D 0 and E Q ( hE T (p)) = E 0 if and only if E Q ( hD T (p)) = D 0 .This condition is called market pricing (MP).

  Du du -( LT ) -,

T0

  du du ≥ D 0 a.s. Actually, we have the following Lemma 4.5. The condition hD T (p) ≥ D 0 holds a.s. if and only if A 0 + T 0 θ u d Su -T 0 du du ≥ D 0 a.s. and, under this condition, (MP) holds if and only if p=0. Proof. Suppose that hD T ≥ D 0 . If hD T = A 0 + T 0 θ u d Su -T 0 du du, then the inequality hD T ≥ D 0 implies that A 0 + T 0 θ u d Su -T 0 du du ≥ D 0 . Otherwise, A 0 + T 0 θ u d Su -T 0 du du ≥ hD T = D 0 + T 0 p u Du du ≥ D 0 . The reverse implication is trivial. At last, under this last condition, note that hD T (p) ≥ hD T (0) = D 0 . Therefore, by (MP), E Q ( hD T (p)) = D 0 implies that T 0 p u Du du = 0 a.s., i.e. p u = 0 du-a.e. 2

T 0 pu 2 Proposition 4 . 8 .

 0248 Du du -γ, as the difference between the l.h.s. and the r.h.s. of the equality is non negative with a zero expectation. This implies that γ ≤ Du du, a.s, i.e. LT = ÃT -DT ≤ 0 a.s. Reciprocally, if ÃT ≤ DT , by (4.11) we get that hD T (p) = ÃT + T (p) does not depend on p and satisfies E Q ( hD T (p)) = D 0 . Suppose that E Q T 0 dr dt < E 0 . Then, there exists a unique credit risk premium p * ∈ R + such that (MP) holds. Proof. We use the same notations than in the proof of Proposition 4.7. Let us introduce the function φ : p ∈ R + → E Q ( hD T (p)) -D 0 . We have φ(p) = E Q T 0 p u Du (p)du -T 0 p u Du (p)du -γ + .

Proposition 4 . 9 . 0 p

 490 Suppose that (MP) holds and let p * ∈ R + be such that φ(p * ) = 0. Then, p * = 0 if and only if ÃT ≥ DT a.s. In the following, we propose to implement an estimation of the credit risk premium in the case where k is a deterministic function. Recall that Dt = D 0 e pt-t 0 kudu t ∈ [0, T ]. In particular, T Dt dt = D 0 p T 0 e pt-t 0 kudu dt =: g(p).

  and (α k ) k≥1 is a positive sequence satisfying

	∞ k=1 α k = +∞ and ∞ k=1 (α

k ) 2 < +∞. Theorem 4.10. The sequence (x k ) k≥1 converges a.s. to x * . Proof. Note that H(v, x) = x1 x≤v + v1 x>γ , in particular x → H(γ, x) is non decreasing. Let us define the non decreasing function C

[START_REF] Schwert | Is Borrowing from Banks More Expensive than Borrowing from the Market? Working Paper[END_REF] gives the example of the average recovery rate for loans that is 80% whereas it is 40% for bonds, which implies that loan credit spreads should be smaller than bond credit spreads. From these aggregated statistics, he argues that the Duffie and Singleton (1999) model predicts that bond spreads should be approximately three times as large as loan spreads. More precisely, when the firm is close to default, bond spreads are significantly higher than loan spreads, but when it is far from default, the loan and bond spreads are similar.

See Papin (2013)'s doctoral dissertation for an interesting discussion on that topic.

In term loans, firms receive generally the full loan amount upfront and repay the loan at maturity, 5-8 years after loan origination, while other term loans are amortized until maturity. Credit lines are more frequent and more complex as they provide contingent liquidity, i.e. they do not draw down the committed loan amount, but keep the credit line as an insurance. The pricing structure of credit lines is more complex as it consists of various fees in addition to the loan spread.


Anirudh Asuri Mukundan
email: anirudh.mukundan@coria.fr

T Menard

A Berlemont
email: alain.berlemont@coria.fr

Jorge César

Brändle De Motta

Brändle De

T Ménard

J C Brändle De Motta

Numerical study of interface reconstruction method in under-resolved regions of the flow for liquid jet primary breakup

Keywords: multiphase, moment of fluid, center of mass, volume fraction

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Atomization of injected liquid fuel in gas turbines and internal combustion (IC) engines play a major role in the production of harmful pollutant emissions. This is due to the fact that the efficiency of liquid fuel combustion and production of pollutants has a direct dependency on the mixing between fuel and oxidizer, which in itself results from the cascade of the mechanisms initiated by the process of atomization of injected liquid fuel stream. Due to the multiphase characteristic and multi-physical aspect of the liquid fuel atomization, and often due to the presence of highly turbulent environment, experimental investigations proves to be a challenging endeavor. This motivates the need for the development of numerical methods to study the atomization process.

An obvious requirement of these methods is least numerical error which will make it more trustable and suitable to be used in detailed simulations such as direct numerical simulations (DNS) of liquid atomization [1][2][3]. In addition, these numerical methods ideally must require less computational resources, thus not taxing even the most powerful supercomputers available today [4].

In the presence of complex topological structures encountered in liquid fuel atomization, the underlying equations are stiff, thus the development and improvement of numerical methods to treat such structures has been an active area of research. Of the many sharp interface reconstruction methods, volume of fluid (VOF) method [5][6][7][8][9], level set (LS) method [10][11][12][13][14], and their variant such as coupled level set volume of fluid (CLSVOF) method [1,3,15] have been proved to perform well on variety of canonical test cases. Although these validated methods prove to be successful in multiple applications, the geometrical property computations can still be inaccurate especially in the under-resolved regions of the flow. This results from the inaccurate computation of unit normal of the interface in such regions. To mitigate this inaccuracy, a numerical interface reconstruction method called moment of fluid (MOF) method [16] is developed in the in-house code ARCHER [1,10,17] and presented in this work. MOF method has been proved -to preserve the orientation of the liquid/gas interface in the context of multiphase flows [16,18]; to be second-order accurate in spatial resolution [19].

The objectives of this paper are: to present MOF method for interface reconstruction in the context of incompressible multiphase flows as a potential interface tracking technique, compare and contrast the results of interface reconstruction for various academic test cases between MOF and CLSVOF methods, and to emphasize the relevance of the usage of the MOF method especially in the under-resolved regions flow regions.

Moment of Fluid (MOF) method

MOF method is an extension of VOF method for interface tracking in the context of multiphase flows. VOF method uses only liquid volume fraction F (0 th moment of the liquid volume) for tracking the interface in every * Corresponding author: anirudh.mukundan@coria.fr 1 mixed computational cell. As a step further, MOF method tracks both the liquid volume fraction and coordinates of the liquid center of mass (COM) x COM (1 st moment of liquid volume) for the interface reconstruction in each such cell. A mixed computational cell, within this study, is defined as the cell in which 0 < F < 1 holds. The definitions of the 0 th and 1 st moments of the liquid volume are given as

F = ω dx Ω dx , (1)
x COM = ω xdx Ω dx , ( 2
)
where ω is the domain of the liquid packet inside the computational cell Ω. The availability of these two parameters establishes a self-sufficiency of the required information to reconstruct the approximate interface in a cell thus, eradicating data requirement from its neigbours. This results in a uniform treatment of the internal and boundary cells in the mesh thus, yielding the resolution of the interface as high as that of the computational mesh itself. Within this work, the MOF method is developed in the in-house code ARCHER.

Interface Reconstruction

The commonly used piecewise linear interface calculation (PLIC) is employed in this work for approximating the original interface. Thus, the equation of the approximated interface in 3D is given by ax + by + cz + d = 0, which represents the equation of a plane. In 2D, this equation becomes the equation of a line. The unit normal n of the interface is therefore expressed as n = [a, b, c] T . The shortest distance of the interface from the cell center is characterized by the parameter d.

The interface reconstruction in MOF method is a constrained optimisation problem as described in [18] wherein n and d have to be determined simultaneously by the satisfaction of the two conditions

| F ref -F act (n, d) |= 0, and (3)
E MOF (n, d) = min
Eq. ( 3) holds

x ref COM -x act COM (n, d) 2 . (4)
All the variables containing the superscript "ref" represents the variables pertaining to the original (reference) interface while those containing the superscript "act" represents the variables pertaining to the reconstructed (actual) interface. The explanation of these conditions is assisted through Figure 1. The shortest distance of the interface 3). To this end, the linear approximated interface (shown by dashed line in Figure 1) is constructed such that the volume of liquid is exactly the same under reference (shown by solid line in Figure 1) and actual interface upto the machine precision. The computation of optimal value of d is carried out using Newton-Raphson iterative method as described in Ménard et al [1]. The interface unit normal is then computed by minimising the error E MOF (also called distance defect) between the coordinates of the reference and the acutal COM of the liquid in the computational cell (c.f. Equation ( 4)).

In order to solve this minimisation problem, n is parameterised [18] using the polar coordinates as follows

n =         a b c         =         sin Φ cos Θ sin Φ sin Θ cos Φ         . (5)
Thus, the Equation ( 4) transforms into a non-linear least square problem for (Φ, Θ), i.e., finding (Φ * , Θ * ) such that the error E MOF is minimum, i.e.,

E MOF (n, d) = g(Φ * , Θ * , d) 2 = min Eq.(3) holds g(Φ, Θ, d) 2 , (6)
where g(Φ, Θ, d) = x ref COM -x act COM . Equation ( 6) is solved numerically for (Φ * , Θ * ) using Gauss-Newton method. Once the optimal parametric angles (Φ * , Θ * ) are known, the components of the unit normal [a, b, c] T can be retrieved using Equation (5). The new reference COM in each cell takes the value of the acutal COM pertaining to the reconstructed interface using which the new values of the normal components are computed.

Interface Advection

The advection of interface in the case of MOF method involves advection of both reference liquid volume fraction F ref and coordinates of reference COM x ref COM . A directionally split numerical scheme is employed for the advection of both quantities. The transport equation and numerical scheme implemented for the advection are presented in the following subsections. Since the advection pertains solely to the reference quantities, the superscript "ref" will be dropped in volume fraction and the coordinates of the COM hereon.

Volume Fraction

The advection equation for the volume fraction solved in ARCHER is given as

∂F ∂t + ∇ • (F u) = c (∇ • u); c = 1, F > 0.5 0, otherwise (7)
The algorithm proposed in the work of Weymouth and Yue [20] is implemented for the advection of the volume fraction. For more details on the finite difference discretization of the above equation and computation of volume fluxes across cell faces, the reader is referred to [20].

Centre of Mass

The advection equation of the reference COM is given as

d dt (x COM ) = u(x COM ), (8)
where u(x COM ) is the velocity field interpolated linearly from the face-centers of the computational cell to the location of the COM. This velocity is non-dimensionalized using the mesh spacing and the time step size ∆t = t n+1 -t n , hence the velocity component becomes the local Courant-Friedrichs-Lewy (CFL) number. In this study, the COM is considered as a Lagrangian particle [16, Appendix A] that is associated with a liquid packet/parcel. An Eulerian Implicit-Lagrangian Explicit (EI-LE) scheme [START_REF] Bnà | 11th International Conference of Numerical Analysis and Applied Mathematics[END_REF] is employed in this study for solving Equation (8). A first order integration scheme of this equation keeping constant velocity over the time step size ∆t yields,

x n+1 COM = x n COM + u(x * COM ) (9)
in which the mode of the scheme is Eulerian Implicit if x * COM = x n+1 COM and Lagrangian Explicit if x * COM = x n COM . In order to have the consistency between the advection of liquid volume fraction and COM, the mode of the scheme for COM advection is switched between EI and LE for a Cartesian direction at each time step, i.e., if x-direction advection of COM for t n → t n+1 is carried out using EI mode, then the x-direction advection of COM for t n+1 → t n+2 will be carried out in LE mode.

In this study, the coupling between the advections of the volume fraction and that of COM is carried out similar to that in the study of Jemison et al [18].

Results and Discussion

In this section, the accuracy, order of convergence, and utility of the MOF interface reconstruction method are analyzed using different canonical test cases. The choices of the validation tests considered in this study are made so that the errors due to reconstruction and advection algorithms used in this study can be assessed. In each of the presented tests, the reconstruction accuracy of MOF and CLSVOF algorithms are compared along with the order of convergence of the error with respect to the spatial resolution of the mesh.

Since the advection errors are almost always accompanied with the reconstruction errors, these errors are reported and addressed using separate test cases for better evaluation of the performance of the reconstruction and advection algorithm implemented in this study. In all the test cases presented in this work, a constant CFL number of 0.5, periodic boundary conditions along xand y-directions, and ∆x = ∆y = const. are used. The time step size ∆t is computed from the CFL number and the mesh spacing.

Translation test

Following the work of Harvie and Fletcher [START_REF] Harvie | [END_REF], the advection of circular liquid droplet of radius r = 0.25 units placed at the center of a 1 × 1 domain is considered for analysing the errors due to reconstruction from MOF and CLSVOF methods. The CLSVOF method developed by Ménard et al [1] is used in this study for the comparison of results with MOF method. The advection is performed along x-direction with unit velocity such that the interface crosses the boundary and comes back to the same initial location at the final time step. Since the velocity is uniform, obeys continuity equation and does not have any gradient in the domain, the errors arising at the end of advection is purely from the reconstruction algorithm. This error is computed as

E = i,j A ij | F nT i,j -F 0 i,j |, ( 10
)
where A i,j is the area of the cell, F 0 and F nT are the liquid volume fractions at the initial and final time iteration instants respectively. The mesh resolutions considered in this test are 16 × 16, 32 × 32, 64 × 64, 128 × 128, and 256 × 256. blue line: CLSVOF in middle subfigure) of the circular liquid droplet after the advection for 32 × 32 mesh. The circle comes back to its initial location after 128 time advancement iterations for this mesh resolution. For such interface topology, it can be seen qualitatively that both MOF and CLSVOF methods have performed well. From the error convergence plot shown in Figure 2c, it is quite obvious that both MOF and CLSVOF display secondorder convergence rate.

Rotation Test

The same case as described in Rider and Kothe [23] is considered in which a circle of radius r = 0.15 units is placed in a 1 × 1 domain with center at (0.5, 0.75), rotates around the center of the domain in a spatially varying field given by u = π 3.14 (0.5 -y), and

v = π 3.14 (x -0.5). (11)
In the exact solution, the circle retains its shape thus, this test assesses the efficiency of the EI-LE advection algorithm implemented in this study. The error in the numerical solution is computed using Equation (10) where nT is the final time iteration at which the interface completes one full revolution around the center of the domain.

The results of the error for the mesh resolutions 32 × 32, 64 × 64, and 128 × 128 are presented and compared in the Figure 3 for MOF and CLSVOF reconstruction algorithms. The initial interface is depicted by black line, the final interfaces from MOF and CLSVOF reconstruction methods are depicted by red and blue lines respectively. It can be observed from this figure that for such a simple advection problem, the advection schemes employed for the volume fraction and COM in this study for MOF method is qualitatively more accurate than CLSVOF method. Furthermore, from Figure 3d for mesh resolutions from 16 × 16 to 128 × 128, it is observed that MOF method exhibits second-order convergence rate and produces lower error than CLSVOF method.

Zalesak's notched disk test

In this test, a slotted circle of fluid rotates around the center of the domain in a solenoidal velocity field. This test is indeed an assessment of the efficiency of the reconstruction algorithm than that of the advection algorithm.

Within this test, a slotted circle of radius r = 0.15 units having slot width of 0.06 units and slot length of 0.2 units is placed in a 1 × 1 domain with center of the circle located at (0.5, 0.75). The rotation of this slotted circle is accomplished using velocity fields given in the Equations ( 11) and (12). Since this velocity is divergencefree throughout the domain, the interface of the Zalesak's notched circular disk is expected to retain the original shape at the end of rotation. The test ends when the first (anticlockwise) revolution of the slotted circular disk is completed. The resulting error between the exact and the reconstructed interface is then computed as

E = i,j | F nT i,j -F 0 i,j | i,j | F 0 i,j | , ( 13
)
where the definitions of F 0 and F nT remains the same as explained in the previous subsection. The results from this test are shown in Figure 4 for 32 × 32, 64 × 64, and 128 × 128 mesh resolutions. The initial interface is depicted using black solid line, the final interface from MOF algorithm is depicted using red solid line while that from CLSVOF algorithm is depicted using blue solid line. 4c, that the error in reconstruction is mainly concentrated in the regions of high curvature of the interface. Furthermore, it is obvious from from Figure 4d that the MOF method is

14 th ICLASS 2018
Interface Reconstruction using MOF method: Numerical Study quantitatively more accurate in capturing the interface under a solenoidal velocity field than CLSVOF method.

Vortex in a box test

This test assesses the ability of an interface reconstruction method to represent the under-resolved structures such as ligaments in a robust manner. In this test, a circular disk is made to undergo high amount of deformation under a given velocity field. The following time reversing velocity field is prescribed over the whole 1 × 1 domain u = -2 sin 2 (πx) sin(πy) cos(πy) cos(πt/T )

v = 2 sin 2 (πy) sin(πx) cos(πx) cos(πt/T ) ( (14)
) 15
where T is the time period of reversal of the velocity field. In this study, T = 6 is chosen for the analysis of the interface reconstruction accuracy. This velocity field stretches and tears the initially circular fluid body as it becomes progressively entrained by the vortex and comes back to its original shape at time t = T . The entrainment is demonstrated as long thin fluid filament spiraling inward towards the vortex center. Unlike Zalesak's notched disk test, the velocity field in this test is non-linear deeming this test to be a more realistic assessment for the reconstruction method.

The results of the interfaces are shown in Figure 5 for the mesh resolutions 64 × 64, 128 × 128, and 256 × 256 with the reference solution (depicted by black solid line in the subfigures) obtained on a 1024 × 1024 grid. The top row in this figure pertain to the time instant t = T /2 corresponding to the maximum stretching of the fluid body and those at the bottom row pertain to the time instant t = T corresponding to the return to the original shape. The reconstruction algorithm proves to be robust in the smooth regions of the interface and it can be said qualitatively that the results become better with the increase in mesh resolution at both the time instants of maximum stretch and at final time instant at which the interface returns back to its original shape. Figure 5d quantifies the reconstruction error computed at time t = T when the interface returns to its initial shape. This error is computed according to Equation (10). From this figure, it can be seen that both MOF and CLSVOF reconstruction algorithms are of the same level of accuracy in terms of reconstruction error. Once again, it can be seen that the MOF method exhibit second-order accuracy in spatial resolution.

Deformation of circular fluid body test

A complex time reversing velocity field [24] given by u = sin(4π(x + 1/2)) sin(4π(y + 1/2)) cos(πt/T ) (16) v = cos(4π(x + 1/2)) cos(4π(y + 1/2)) cos(πt/T ) (17) is used in this test for the deformation of a circular bubble under extreme conditions. This velocity field will induce radical deformation in the topology of the fluid body thereby providing a more stringent test for the reconstruction method than those presented before. To this end, a circular fluid body of radius r = 0.15 units is placed with its center at (0.5, 0.75) in a 1 × 1 domain. At its maximum deformation, the circular fluid body is entrained into two of the four nearest vortices with a small portion (thin filaments) entrapped by the two of the nearby vortices. The time period for velocity reversal is chosen as T = 2.0 for this test. The results are shown in Figures 6b to 6d for the 128 × 128, 256 × 256, and 512 × 512 mesh resolutions for t = T /2 (top row) and t = T (bottom row) with initial interface (black line), reconstructed interfaces using MOF (red line), and CLSVOF (blue line) for the corresponding time instants. In addition, a zoomed view of the miniature ligament is found in the final time iteration of the 512 × 512 mesh is shown. From this zoomed view along with the underlying mesh, it can be seen that the depicted interface structure is extremely under-resolved. Thus, even with MOF method spurious oscillations in the final interface are seen in the figures in the top row of Figures 6b to 6d. Albeit this current limitation, even under extreme velocity fields which is tearing the interface apart, the MOF is observed to perform satisfactorily for interface reconstruction under extreme deformations of the liquid body topology (c.f. top row in Figures 6b to 6d).

The error from the reconstruction method is quantified using the Equation (10). The plot of the convergence of this error for 32 × 32, 64 × 64, 128 × 128, 256 × 256, and 512 × 512 mesh resolutions is shown in Figure 6a. It can be seen that the MOF method is performing to a satisfactory extent for any mesh resolution considered in this test case considering the fact that the interface is entirely deformed due to the prescribed velocity field. Furthemore, it can be seen that both the MOF and CLSVOF method exhibits a third-order convergence rate for all the spatial resolutions considered for this test case. The reason for this behavior in comparison to the other test cases could be attributed to the fact that as the mesh resolution increases, particularly for this test case, smaller structures of the interface are increasingly identified thereby enhancing the convergence rate of the error.

Conclusions

Study of liquid fuel atomization process is essential in order to meet the stringent emission norms. This is due to the fact that this process has a direct influence on the pollutant emissions. Numerical study of atomization under turbulent conditions gives insights which are otherwise challenging to obtain through experimental investigatations.

Within this work, MOF method was presented for reconstruction of liquid/gas interface in the context of incompressible multiphase flows to the application of numerical analysis of primary breakup of liquid fuel. This method involves computation of interface unit normal by reducing the defect in the first order moment of liquid volume in a volume conservative manner. For various test cases presented in this study, MOF method proved to be more accurate than and at least as accurate as that of CLSVOF method in terms of reconstruction error and preservation of the shape and orientation of the interface topologies. Furthermore, MOF method demonstrated its superior capability of preserving the interface for topologies found in under-resolved flow structures such as thin filaments and ligaments. This deems MOF method as an interface reconstruction method to be used in DNS for studying primary atomization.

In future, development of a hybrid MOF-CLSVOF method for interface reconstruction is a promising improvement thereby combining the advantages of both methods for interface reconstruction. MOF method requires more computational resources than CLSVOF method for the interface reconstruction mainly due to its unit normal computation from the COM. Thus, it is more prudent to use CLSVOF for interface reconstruction especially for the topologies that are far less complicated making it more sensible to develop a hybrid MOF-CLSVOF method. Also, the volume fraction and COM of gas phase along with those pertaining to the liquid phase should be considered for optimal interface reconstruction in the next studies.

Figure 1 :

 1 Figure 1: Mathematical Representation of MOF method

Figure 2 :

 2 Figure 2: 32 × 32 mesh; Black line (nt = 0): Initial; Red line (nt = 128): MOF; Blue line (nt = 128): CLSVOF

Figure 3 :

 3 Figure 3: Initial and final interfaces of circle after one full rotation and error convergence plot

Figure 4 :

 4 Figure 4: Initial and final interface of Zalesak's slotted disk after one full revolution and error convergence plot

  Error convergence for vortex in a box test

Figure 5 :

 5 Figure 5: Vortex in a box test; t = T /2 (top row) and t = T (bottom row) and error convergence plot

10 -3 10 -

 10 Error convergence for deformation of circular fluid body (b) 128 × 128 (c) 256 × 256 (d) 512 × 512

Figure 6 :

 6 Figure 6: Interface reconstruction for deformation of circular fluid body t = T /2 (top row) and t = T (bottom row) and error convergence plot

Acknowledgements

The funding for this project from the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No. 675676 is gratefully acknowledged.


M Usman
email: musmanch@yahoo.com

K Hanna

P Faure

Science of The Total Environment Remediation of oil-contaminated harbor sediments by chemical oxidation

Keywords: Oil hydrocarbons, Harbor sediments, Remediation, Fenton-like, Persulfate, Magnetite

Remediation of oil-contaminated harbor sediments by chemical oxidation

Introduction

Oil hydrocarbons are widespread pollutants in sub-surface sediments which pose serious threats to terrestrial and aquatic environment. Therefore, development of novel and effective remediation techniques for oil-contaminated sediments is of considerable environmental interest [START_REF] Lim | A comprehensive guide of remediation technologies for oil contaminated soil -Present works and future directions[END_REF]. Various remediation techniques have been developed including biological [START_REF] Chaîneau | Bioremediation of a crude oil-polluted soil: Biodegradation, leaching and toxicity assessments[END_REF][START_REF] Benyahia | Bioremediation of crude oil contaminated soils: A black art or an engineering challenge? Process Saf[END_REF], chemical (Usman et al., 2012a;2013;[START_REF] Chen | Remediation of diesel-contaminated soil using in situ chemical oxidation (ISCO) and the effects of common oxidants on the indigenous microbial community: a comparison study[END_REF] and combined approaches [START_REF] Silva-Castro | Bioremediation of diesel-polluted soil using biostimulation as post-treatment after oxidation with Fenton-like reagents: Assays in a pilot plant[END_REF][START_REF] Souza | Biosurfactant-enhanced hydrocarbon bioremediation: An overview[END_REF][START_REF] Agarwal | Remediation technologies for oil-contaminated sediments[END_REF]Xu et al., 2016a). Bioremediation is the most commonly used technique to remediate oil contaminated matrices (soils, sediments and water) primarily because of its cost-effective and environment friendly nature [START_REF] Lim | A comprehensive guide of remediation technologies for oil contaminated soil -Present works and future directions[END_REF]. Majority of molecules in crude oil spills and refined products are biodegradable and thus can largely be removed by bioremediation [START_REF] Prince | Petroleum spill bioremediation in marine environments[END_REF]. However, bioremediation often leaves nonbiodegradable residues due to incomplete mineralization and thus cannot meet the stringent environmental cleanup standards reported by environmental regulations [START_REF] Płaza | The application of bioassays as indicators of petroleum-contaminated soil remediation[END_REF].

Contamination level had reached a plateau for long time in soil subjected to bioremediation for two years [START_REF] Lu | Removal of residual contaminants in petroleum-contaminated soil by Fenton-like oxidation[END_REF]. Bioremediation or natural attenuation of n-alkanes present in oil lead to their degradation whereby leaving unresolved complex mixture (UCM) composed of iso-and cyclo-alkanes, which is totally refractory to microbial attack [START_REF] Chaillan | Factors inhibiting bioremediation of soil contaminated with weathered oils and drill cuttings[END_REF][START_REF] Farrington | Unresolved Complex Mixture" (UCM): A brief history of the term and moving beyond it[END_REF]. However, chemical oxidation has shown strong efficiency to remove oil hydrocarbons where different oxidants are used to degrade organic pollutants.

Recent studies employing chemical oxidation for the remediation of oil hydrocarbons are summarized in Table 1. Among chemical oxidation-based treatments, Fenton and persulfate oxidation are showing great remediation potential to remove oil hydrocarbons. Both oxidants require ferrous iron Fe(II) to produce stronger radicals through following reactions: Fe(II) + H 2 O 2 → Fe(III) + HO -+ HO Activation by soluble Fe(II) requires an acidic pH (~ 3) to avoid its precipitation which is harmful for soil quality and microbial population [START_REF] Laurent | Oxidation of a PAH polluted soil using modified Fenton reaction in unsaturated condition affects biological and physicochemical properties[END_REF]. This initial acidification can be avoided by using chelating agents or iron minerals to perform chemical oxidation at circumneutral pH (Usman et al., 2016b). Iron oxides especially mixed-valent Fe(II)-Fe(III) minerals like magnetite [Fe(II)Fe(III) 2 O 4 ] and green rust were found more effective to catalyze chemical oxidation as compared to the ferric minerals or soluble Fe(II) at circumneutral pH due to the presence of structural Fe(II) [START_REF] Usman | Magnetite and green rust: synthesis, properties, and environmental applications of mixed-valent iron minerals[END_REF]. Previous batch (Usman et al., 2012a) and column experiments [START_REF] Usman | Treatment of hydrocarbon contamination under flow through conditions by using magnetite catalyzed chemical oxidation[END_REF] showed that magnetite-catalyzed chemical oxidation (Fenton-like and persulfate) is highly efficient to remove oil hydrocarbons (90% degradation) in crude oil as well its refractory product after biodegradation (UCM). However, these findings were based on artificially contaminated sand where oil hydrocarbons were spiked on sand after extraction from contaminated sediments.

To date, no study has been published on use of chemical oxidation to remediate real oilcontaminated sediments. It should be noted that chemical oxidation resulted in significant removal (80%) of polycyclic aromatic hydrocarbons (PAHs) when extracted from aged contaminated soils and spiked in sand. However, these treatments were ineffective (<5% degradation) in aged contaminated soils (without any extraction) which has been attributed to the contaminant unavailability (Usman et al., 2012b;2012c;2016a;[START_REF] Biache | Effect of preheating on the chemical oxidation efficiency:Implications for the PAH availability measurement in contaminated soils[END_REF]. In Dichloromethane and chloroform, both purchased from VWR were used as received.

Deionized water was produced with a Milli-Q system from Millipore.

Magnetite used in the present study was synthesized and characterized in the context of previous study (Usman et al., 2012a). This magnetite is characterized as nano-magnetite (particle size range of 30 -50 nm) with BET surface area of 103 ± 2 m 2 g -1 , PZC of 7.9 and Fe(II)/Fe(III) ratio of 0.46.

Sediment samples

This study concerned three contaminated sediments from two different harbors, and one petroleum industrial channel. The locations of these three sites are shown in the Figure 1 petroleum industrial channel located in the south of France. These samples were freeze-dried, sieved at 2 mm, placed in glass bottles, and stored at 5 °C before analysis and/or treatment.

Properties of the tested sediments are reported in Table 2.

These three locations are affected by petroleum by-product contamination for several decades and polluted sediments have already undergone natural biological alteration.

Oxidation procedures

Following experiments were performed to evaluate the chemical oxidation of contaminated [START_REF] Liang | Persulfate oxidation for in situ remediation of TCE. II. Activated by chelated ferrous ion[END_REF]Xue et al., 2009b;Usman et al., 2012a;2012b), we used higher dose of H 2 O 2 than persulfate. Blank experiments were performed by using magnetite alone without any oxidant but no pollutant degradation was observed (data not shown). Magnetite was added in concerned sediments to achieve 10% w/w as final contents (0.1 g of magnetite for 1 g of sediment). Equivalent molar amount of Fe was used to compare the efficiency of both catalysts (soluble Fe(II) and magnetite) to activate the tested oxidants. Standard procedure is detailed in our previous study by Usman et al. (2012a). Sediments (2 g) were stirred for 15 minutes with suitable amount of water to have a final volume of 20 mL after oxidant addition. Photocatalytic degradation was avoided by using aluminum foil coverings. Theses slurries were reacted for one week under vigorous magnetic-stirring. Reaction was then stopped by freezing these batches. After two days, water was removed by freeze drying these samples.

One week time point was chosen based on previous studies where no further degradation was observed after one week (Usman et al., 2012a;2012b) while oxidants were fully consumed after this period as confirmed by iodometric titration [START_REF] Rybnikova | Removal of PCBs in contaminated soils by means of chemical reduction and advanced oxidation processes[END_REF]. All experiments were performed in duplicates and results are expressed as a mean value of the 2 experiments and standard deviation of the duplicates was less than 5%.

Extraction and analysis

Quantification of extractable organic matter (EOM): Freeze-dried sediments were extracted with chloroform at 60 °C for 45 min as detailed previously (Usman et al., 2016a). This was performed in a glass vessel which was connected to a water-cooled reflux condenser to let reflux the evaporated solvent. Temperature of extraction medium was homogenized by magnetic stirring. After extraction, the vessel was allowed to cool at room temperature before opening it. The supernatant was filtered through a pre-washed glass microfiber filter (d = 45 mm, porosity = 0.7 mm; Whatman). Volume of organic extract was reduced to 20 mL under N 2 and 5 mL of this solution was transferred in a pre-weighed vial. The amount of extractable organic matter (EOM) was determined by weighing the vial after evaporation of the solvent.

This method resulted in complete recovery of EOM in previous studies (Usman et al., 2012a;2012c).

Hydrocarbon Index (HI) quantification:

The HI was measured according to the ISO 16703:2004 procedure using a GC-FID 7890 Agilent technologies. Briefly, a standard mixture (ASTM D5307) is used to define the integration limits of the chromatogram area.

These limits correspond to the n-decane (C 10 H 22 ) and n-tetracontane (C 40 H 82 ) retention time.

An ISO 11046 standardized oil is used as standard solution (external calibration). Seven solutions of different concentrations (0.5, 1, 2, 4, 6, 8 and 10 mg mL -1 ) are prepared from this standardized oil solution and injected into the chromatograph in order to obtain calibration curves. The integration of the area under each chromatogram is corrected by the background noise (by integrating the area chromatogram obtained with a sample of pure dichloromethane that serves as a control).

Results and discussion

Characterization of sampled sediments

The main properties of tested sediments are reported in Observed order of TOC contents could be correlated to increasing silt contents and decreasing sand contents as reported by [START_REF] Augustin | Relationships between soil carbon and soil texture in the Northern Great plains[END_REF]. Indeed, total silt was positively correlated (P ≤ 0.01) to the TOC while sand was negatively correlated (P ≤ 0.10) to it. An order similar to TOC was observed for HI values: LGH (3.7 g kg -1 ) < TMH (5.5 g kg -1 ) < PIC (9.0 g kg -1 ). These sediments contained higher CaCO 3 amounts of 18.63% (LGH), 28.43 (TMH) and 23.1 (PIC). Higher content of Cu, Pb and Zn were also detected in these sediments (Table 2).

The molecular distributions of initial sediments (Figure 2) are dominated by UCM which comprise a mixture of iso-and cyclo-alkanes [START_REF] Tolosa | Aliphatic and aromatic hydrocarbons in coastal caspian Sea sediments[END_REF] and cannot be resolved by classical chromatographic column. The main identified products correspond to pentacyclic triterpanes. Such signatures are typical of biodegraded and weathered oil containing compounds which are resistant to biodegradation [START_REF] Tolosa | Aliphatic and aromatic hydrocarbons in coastal caspian Sea sediments[END_REF][START_REF] Farrington | Unresolved Complex Mixture" (UCM): A brief history of the term and moving beyond it[END_REF]. Sampled areas were under natural attenuation for years which eliminated the easily biodegradable n-alkanes leaving UCM as the sole contamination in three sediments. GC-FID chromatograms (Fig. 2) provide a visual representation of efficiency of chemical oxidation to degrade UCM as observed previously in spiked sand (Usman et al., 2012a).

LGH Evolution of extractable organic matter (EOM) was obtained by measuring the weight of organic extract before and after chemical oxidation (Figure 3). Obtained results (Figure 3) indicated a significant decrease in EOM (25 -50%) when magnetite was used as a catalyst. In TMH, EOM removal extent was lowest (~25% by both oxidants) that might be correlated to its higher carbonate contents (28.43% in TMH sediments as compare to the 18.63% and 23.1% in LGH and PIC sediments, respectively).

Chemical oxidation of contaminated sediments

Carbonates were shown to have scavenging effects for radicals and thus resulting in decreased oxidation efficiency [START_REF] Liang | Influences of carbonate and chloride ions on persulfate oxidation of trichloroethylene at 20 °C[END_REF][START_REF] Grebel | Effect of halide ions and carbonates on organic contaminant degradation by hydroxyl radical-based advanced oxidation processes in saline waters[END_REF]. Loss of EOM was 50% and 32% in LGH sediments and 48% and 35% in PIC sediment by H 2 O 2 /magnetite and Na 2 S 2 O 8 /magnetite, respectively. However, use of both oxidants alone or with soluble Fe(II) resulted in lower removal of EOM (< 5 %, data not shown) that can be associated to the inability of Fe(II) to catalyze chemical oxidation at circumneutral pH (Xue et al., 2009b;Usman et al., 2012a). Soluble Fe(II) precipitates at circumneutral pH which obstructed its catalytic ability. Reaction pH is very important factor controlling the remediation efficiency.

Owing to its strong stability, magnetite was found effective at a wide range of pH (Usman et al., 2016b).

Hydrocarbon index (HI) measured by GC-FID is presented in Figure 4. [START_REF] Xu | Effect of Fenton pre-oxidation on mobilization of nutrients and efficient subsequent bioremediation of crude oilcontaminated soil[END_REF]. It should, however, be noted that total concentration of native Fe(II) in their soil was 4.8 g kg -1 which could explain the reported oxidation efficiency without any added Fe(II) by [START_REF] Xu | Effect of Fenton pre-oxidation on mobilization of nutrients and efficient subsequent bioremediation of crude oilcontaminated soil[END_REF]. Native Fe in contaminated soils has been found to catalyze Fenton-like oxidation with no Fe addition for the removal of PAHs in contaminated soils with 57.5 g kg -1 endogenous Fe oxides (amorphous + crystalline) [START_REF] Watts | Mineralization of a sorbed polycyclic aromatic hydrocarbon in two soils using catalyzed hydrogen peroxide[END_REF]. Similary, [START_REF] Flotron | Removal of sorbed polycyclic aromatic hydrocarbons from soil, sludge and sediment samples using the Fenton's reagent process[END_REF] bearing 94 and 30 g kg -1 as total native Fe, respectively. Owing to the presence of native Fe (16.4 g kg -1 ), 52% of PAHs were degraded by H 2 O 2 (0.4 g g -1 soil) under saturated conditions without any Fe(II) addition [START_REF] Palmroth | Effect of modified Fenton's reaction on microbial activity and removal of PAHs in creosote oil contaminated soil[END_REF]. Role of endogenous Fe (18.2 g kg -1 ) was also highlighted by [START_REF] Pardo | Optimization of the application of the Fenton chemistry for the remediation of a contaminated soil with polycyclic aromatic hydrocarbons[END_REF] when they observed similar degradation of PAHs by Fenton oxidation (882 mM H 2 O 2 without pH adjustment) by adding 1 mM Fe(III) or without external addition of Fe. However, Usman et al. (2012b;2012c) reported that native iron (25 g kg -1 Fe 2 O 3 ) was unable to catalyze H 2 O 2 and Na 2 S 2 O 8 to remediate PAH contaminated soils at circumneutral pH. Similarly, [START_REF] Ahmad | Persulfate activation by subsurface minerals[END_REF] observed that native iron oxides (8.8 g kg -1 amorphous + crystalline) in soil did not promote measureable activation of Na 2 S 2 O 8 in natural soils. They suggested that quantities of native iron oxides are not sufficient to catalyze persulfate oxidation in soils. However, even at varying concentrations of native iron (15 -110 g kg -1 ) in tested sediments, it was found ineffective in promoting the oxidation reaction for harbor sediments. Thus, catalytic role of native Fe oxides is variable depending on many factors including their nature, total content and availability of Fe oxides. A good correlation was observed between hydrocarbon index and extractable organic matter for H 2 O 2 /magnetite highlighting the non-specific nature of chemical oxidation to degrade the non-target compounds. However, Na 2 S 2 O 8 /magnetite resulted in lower removal of EOM while higher loss of HI especially in the PIC sediments. It should be noted that HI quantifies only low molecular weight compounds (less than 40 carbons), whereas EOM quantifies low as well as high molecular weight compounds. This observed difference in degradation of EOM and HI by Na 2 S 2 O 8 /magnetite could be caused by selective behavior of persulfate oxidation as observed previously where it showed less efficiency towards high molecular weight PAHs while Fenton oxidation was equally efficient for all PAHs (Usman et al., 2012b;2012c). Similarly, when crude oil rich in n-alkanes with pristane and phytane spiked on sand was subjected to magnetite-catalyzed chemical oxidation (Usman et al., 2012a), H 2 O 2 /magnetite showed non-selective oxidation with equal efficiency towards all n-alkanes (C 21 -C 36 ). On the other hand, Na 2 S 2 O 8 /magnetite was less reactive towards high molecular weight n-alkanes (C 31 -C 36 ). However, it is difficult to point out such selectivity in studied sediments because they are contaminated with UCM whereby individual components of this complex mixture are difficult to resolve [START_REF] Farrington | Unresolved Complex Mixture" (UCM): A brief history of the term and moving beyond it[END_REF]. Observed difference in term of degradation of EOM and HI for persulfate for PIC sediments could also be correlated to the lower affinity of Na 2 S 2 O 8 towards natural organic compounds leading to higher removal of target pollutants [START_REF] Lim | A comprehensive guide of remediation technologies for oil contaminated soil -Present works and future directions[END_REF]. Very recently, non-productive consumption of Na 2 S 2 O 8 was lower than 10% after 56 days of reaction in non-polluted soil (despite 6.5% CaCO 3 ) highlighting the low reactivity of oxidant with different soil components [START_REF] Lominchar | Remediation of aged diesel contaminated soil by alkaline activated persulfate[END_REF]. Rather, persulfate was also found to be activated at basic pH by organic compounds similar to those present in soil organic matter [START_REF] Ahmad | Persulfate activation by subsurface minerals[END_REF]. Similarly, the impact of soil matrix was less pronounced on the efficiency of persulfate oxidation to remove PAHs from historically contaminated soil as compared to the Fenton oxidation system (Usman et al., 2012b;2012c).

Magnetite-catalyzed chemical oxidation vs. common oxidation systems

As mentioned above, the available literature severely lacks investigations on the chemical oxidation of oil hydrocarbons in real contaminated sediments (Table 1) and thus a fair comparison of present findings with the literature is not straightforward. Previous studies were either performed on artificially contaminated sand/soil or soil obtained from contaminated site (Table 1) and present study is the first report of using magnetite-catalyzed chemical oxidation in real sediments.

Most of the remediation studies to degrade oil hydrocarbons employing H 2 O 2 and Na 2 S 2 O 8 were performed at acidic pH (≤ 3) by Fe-catalyzed chemical oxidation [START_REF] Kong | Treatment of petroleum-contaminated soils using iron mineral catalyzed hydrogen peroxide[END_REF][START_REF] Goi | Ozonation and Fenton treatment for remediation of riesel fuel contaminated soil[END_REF][START_REF] Do | Effect of metal oxides on the reactivity of persulfate/Fe(II) in the remediation of diesel-contaminated soil and sand[END_REF][START_REF] Villa | Diesel degradation in soil by Fenton process[END_REF][START_REF] Pardo | Remediation of a Biodiesel Blend-Contaminated Soil by Using a Modified Fenton Process[END_REF] (Table 1). Due to detrimental effects of this initial acidification, efforts were made to perform chemical oxidation of oil hydrocarbons at circumneutral pH by using catalytic ability of basic oxygen furnace slag [START_REF] Tsai | Treatment of petroleum-hydrocarbon contaminated soils using hydrogen peroxide oxidation catalyzed by waste basic oxygen furnace slag[END_REF], soluble Fe(III) or Fe(II) with chelating agents (CAs) like EDTA or trisodium citrate or citric acid [START_REF] Lu | Removal of residual contaminants in petroleum-contaminated soil by Fenton-like oxidation[END_REF][START_REF] Pardo | Remediation of a Biodiesel Blend-Contaminated Soil by Using a Modified Fenton Process[END_REF]Xu et al., 2016b;[START_REF] Xu | Effect of Fenton pre-oxidation on mobilization of nutrients and efficient subsequent bioremediation of crude oilcontaminated soil[END_REF], magnetite (Usman et al., 2012a;2013;[START_REF] Jamialahmadi | Remediation of a diesel-contaminated soil using a Fenton-like advanced oxidation process: Optimization by response surface methodology[END_REF] and zerovalent Fe [START_REF] Jamialahmadi | Remediation of a diesel-contaminated soil using a Fenton-like advanced oxidation process: Optimization by response surface methodology[END_REF]. Very high pH (>12) also resulted in strong activation of Na 2 S 2 O 8 for the removal of diesel. Use of chelating agents resulted in effective removal of oil hydrocarbons (Table 1). However, these CAs being organic in nature will act as oxidant scavenger by competing with the target pollutants [START_REF] Pardo | Remediation of a Biodiesel Blend-Contaminated Soil by Using a Modified Fenton Process[END_REF][START_REF] Usman | Comment on "A comprehensive guide of remediation technologies for oil contaminated soil -Present works and future directions[END_REF]. For example, [START_REF] Pardo | Remediation of a Biodiesel Blend-Contaminated Soil by Using a Modified Fenton Process[END_REF] applied Fenton oxidation to treat biodiesel contaminated soil with or without trisodium citrate as CA (for details of experimental conditions, see Table 1).

Obtained oil hydrocarbon removal reached 75% in the absence of CA (at pH <3). Use of CA allowed them to perform chemical oxidation at higher pH (5-6) but degradation efficiency decreased to 37% which was associated to the consumption of oxidant by CA at higher pH.

Use of CA can also be pH sensitive as reported by [START_REF] Lu | Removal of residual contaminants in petroleum-contaminated soil by Fenton-like oxidation[END_REF]. H 2 O 2 catalyzed by Fe(III)

in the presence of a CA (EDTA) was more efficient to degrade petroleum hydrocarbons at pH 7 while oxidation efficiency decreased at pH 6 which was associated to quick degradation of EDTA leading to the precipitation of Fe(III). However as reported in literature, magnetite has shown strong stability and catalytic capacity over a wide range of pH and thus, it offers control over drawbacks associated to initial acidification or CAs. Moreover owing to its ferromagnetic nature, it can be easily recovered after reaction by magnetic separation and can be used for further oxidation cycles without any significant decline in catalytic ability (Xue et al., 2009a;[START_REF] Jia | TCE Degradation in groundwater by chelators-assisted Fenton-like reaction of magnetite: Sand columns demonstration[END_REF]. Other reusability studies indicated a decline in oxidation efficiency of 20% after 8 oxidation cycles for the removal of phenol in aqueous solutions [START_REF] Zhang | Superparamagnetic Fe 3 O 4 nanoparticles as catalysts for the catalytic oxidation of phenolic and aniline compounds[END_REF]. However, reusability and/or fate of magnetite has rarely been investigated in soils or sediments.

Regarding efficiency of magnetite-catalyzed chemical oxidation in present study, degradation of hydrocarbons varied between 40 -70% when tested sediments were subjected to chemical oxidation at circumneutral pH with oxidant:Fe molar ratio equal to 20:1 (H 2 O 2 ) and 2:1 (Na 2 S 2 O 8 ) with 10% w/w magnetite. This would be economically relevant by avoiding initial acidification and/or use of chelating agents.

It is worth noting that efficiency of chemical oxidation is strongly limited by the availability of polycyclic aromatic hydrocarbons (PAHs) in aged contaminated soils and sediments [START_REF] Flotron | Removal of sorbed polycyclic aromatic hydrocarbons from soil, sludge and sediment samples using the Fenton's reagent process[END_REF][START_REF] Choi | Degradation of pyrene in cetylpyridinium chloride-aided soil washing wastewater by pyrite Fenton reaction[END_REF]Usman et al., 2016b). As a matter of fact, H 2 O 2 /magnetite or Na 2 S 2 O 8 /magnetite did not cause any PAH degradation unless tested soils were subjected to chemical (Usman et al., 2012c) or thermal pretreatment [START_REF] Biache | Effect of preheating on the chemical oxidation efficiency:Implications for the PAH availability measurement in contaminated soils[END_REF]Usman et al., 2016a). However, no such pollutant unavailability was observed for oil hydrocarbons in tested sediments. It might be correlated to difference between the nature of both matrices (soils vs. sediments), history of contamination and type of pollutants (aliphatic vs. aromatic) [START_REF] Flotron | Removal of sorbed polycyclic aromatic hydrocarbons from soil, sludge and sediment samples using the Fenton's reagent process[END_REF][START_REF] Trellu | Characteristics of PAH tar oil contaminated soils-Black particles, resins and implications for treatment strategies[END_REF][START_REF] Lominchar | Remediation of aged diesel contaminated soil by alkaline activated persulfate[END_REF]. The specific deposition, the preservation in water can probably explain the difference in term of availability between soils and sediments. Moreover, in our case, the contamination originated from petroleum source (a liquid) compared to coal-tar (source of PAHs) which becomes solid with aging while losing the low molecular compounds which leads to enrichment of heavy compounds in coal-tar particles [START_REF] Trellu | Characteristics of PAH tar oil contaminated soils-Black particles, resins and implications for treatment strategies[END_REF]. Coal-tar can be observed in historically contaminated soil as a small ball-like aggregate with low specific area as compared to the petroleum pollution which exists as a thin film on minerals in sediments [START_REF] Trellu | Characteristics of PAH tar oil contaminated soils-Black particles, resins and implications for treatment strategies[END_REF].

In addition to the in-situ application, magnetite-catalyzed chemical oxidation could offer a quick and cost-effective alternative to conventional remediation strategies (landfilling, stabilization etc.) for dredged sediments [START_REF] Akcil | A review of approaches and techniques used in aquatic contaminated sediments: metal removal and stabilization by chemical and biotechnological processes[END_REF] ensuring safe reuse of dredged sediments.

Conclusion

This is the premier study reporting the use of magnetite-catalyzed chemical oxidation to remove biorefractory oil residues (remaining of natural attenuation) in different historically contaminated sediments. Magnetite showed strong efficiency to catalyze chemical oxidation at circumneutral pH. Oxidation efficiency was, however, limited by the nature and properties of sediments. No issue related to pollutant availability was observed as was the case for PAH aged contaminated soils which might be correlated to the difference of contaminant type (aliphatic hydrocarbons and PAHs) and mediums (soils and sediments). Due to the strong efficiency of magnetite to degrade biorefractory oil residues without any initial acidification (circumneutral pH), this study may have strong implications to design an innovative standalone technology for in-situ and ex-situ remediation of contaminated sediments.

  and are detailed below: LGH Sediments were sampled in February 2007 by dredging in Le Guilvinec harbor (north-western of France), TMH Sediments were dredged in Toulon military harbor (southern France) in April 2008, PIC sediments were sampled in 2011 in a

Figure 1 :

 1 Figure 1: Location of tested sediments in France

Figure 2 :

 2 Figure 2: GC-FID spectra of tested sediments before (initial) and after chemical oxidation by H 2 O 2 /magnetite and Na 2 S 2 O 8 /magnetite. Tested sediments were obtained from Le Guilvinec harbor (LGH), Toulon military harbor (TMH) and petroleum industrial channel (PIC). Each organic extract was analyzed by GC-FID using the same procedure (same initial amount of soil before treatment, same amount of solvent for the extraction of sediment after treatment and same volume injected in the GC-FID).

  Contaminated sediments were subjected to chemical oxidation in batch experiments. Efficiency of various oxidation treatments was tested at circumneutral pH including H 2 O 2 alone, H 2 O 2 /Fe(II), H 2 O 2 /magnetite, Na 2 S 2 O 8 alone, Na 2 S 2 O 8 /Fe(II), and Na 2 S 2 O 8 /magnetite.

Figure 3 :

 3 Figure 3: Extractable organic matter (EOM) before (T0) and after chemical oxidation by H 2 O 2 /magnetite (A) and Na 2 S 2 O 8 /magnetite (B) in contaminated sediments obtained from Le Guilvinec harbor (LGH), Toulon military harbor (TMH) and petroleum industrial channel (PIC). Experimental conditions were sediment = 2 g, volume of solution = 20 mL, oxidant doses were used according to the oxidant:Fe molar ratio equal to 20:1 (H 2 O 2 ) and 2:1 (Na 2 S 2 O 8 ) with 10% w/w magnetite.

Figure 4 :

 4 Figure 4: Hydrocarbon index (HI) before (T0) and after chemical oxidation by H 2 O 2 /magnetite (A) and Na 2 S 2 O 8 /magnetite (B) in contaminated sediments obtained from LeGuilvinec harbor (LGH), Toulon military harbor (TMH) and petroleum industrial channel (PIC). Experimental conditions were sediment = 2 g, volume of solution = 20 mL, oxidant doses were used according to the oxidant:Fe molar ratio equal to 20:1 (H 2 O 2 ) and 2:1 (Na 2 S 2 O 8 ) with 10% w/w magnetite.

Table 1 :

 1 Summary of studies using chemical oxidation to remove oil hydrocarbons in contaminated soils.

	Oxidant Catalyst type	Target pollutant and its initial concentration	Medium	Experimental conditions	Extent / rate of pollutant degradation	Ref.
		Magnetite			pH = 3,
	H 2 O 2	(Commercial) Screened through	Petroleum hydrocarbons 1 g Kg -1	Spiked sand/ Batch	magnetite = 5 wt%, H 2 O 2 = 15 wt%	50% hydrocarbon removal after 8 days	(Kong et al., 1998)
		0.3 mm sieve.
					pH = 3,
	H 2 O 2	Soluble Fe(II)	Diesel 10-12 g Kg -1	Spiked sand and peat soil/ Batch	H 2 O 2 :Fe(II) molar ratio = 10:1 liquid/solid = 1 and 20	>80% removal of hydrocarbons after 72 h in sand	(Goi et al., 2006)
					mL g -1 for sand and peat
	H 2 O 2	Basic oxygen furnace slag	Fuel oil or diesel TPH = 10 g kg -1	Spiked sand and sandy loam/ Batch	pH = 6.7, H 2 O 2 = 15% Slag = 100 g kg -1	TPH removal of 96% and 76% in diesel and fuel oil tests after 40 h	(Tsai and Kao, 2009)
	H 2 O 2	Soluble Fe(II) in the presence or absence of natural soil Fe (hematite mainly) minerals	Diesel 5 g Kg -1	Spiked loam soil/ Batch	pH = 2.5 -3, Fe(II) = 12 mM H 2 O 2 = 0.045 g g -1 soil	28% of diesel was degraded of soluble Fe(II) after 84 h with soil minerals (without Fe(II)) while degradation reached 48% in the presence	(Villa et al., 2010)
	H 2 O 2	Soluble Fe(III) with EDTA	Petroleum hydrocarbons Dichloromethane-extractable organics 14.8 g Kg -1	Soil samples Batch collected from areas adjacent to oil wells/	pH = 7, H 2 O 2 = 2.45 mol Kg -1 soil, H 2 O 2 :Fe(III) molar ratio = 200:1.	>80% removal of hydrocarbons	(Lu et al., 2010)
	Na 2 S 2 O 8 Soluble Fe(II)	Diesel 5 g Kg -1	Spiked sand and soil/ Batch	pH = 3, Oxidant:Fe molar ratio = 100:1	38% removal but inclusion of metal oxides improved the oxidation efficiency	(Do et al., 2010)
	Na 2 S 2 O 8 Soluble Fe(II)	Diesel TPH = 3.5 g Kg -1	Soil from an oil refining plant/Batch	Oxidant = 20% Oxidant:Fe = 1:0.01	Total removal 40% after 60 days, Reaction rate constant = 1.47 × 10 -2 day -1	(Yen et al., 2011)
	H 2 O 2 or Na 2 S 2 O 8	Soluble Fe(II) Or Magnetite (lab synthesized)	Petroleum hydrocarbons 4 g Kg -1 of crude and weathered oil	Spiked sand/ Batch	pH = 6.7, 20 -25 °C, magnetite = 10% w/w, 10:1 (H 2 O 2 ) and 1:1 oxidant:Fe molar ratio =	80-90% of hydrocarbon removal by both oxidants in after one week. <10 % of the presence of magnetite	(Usman et al., 2012a)

the present study, we aim to evaluate the efficiency of hydrogen peroxide (H 2 O 2 ) and persulfate (Na 2 S 2 O 8 ) catalyzed by magnetite to remove oil hydrocarbons in aged contaminated sediments. Remediation of three kinds of sediments sampled from different locations in France were investigated. These sediments were subjected to magnetite-catalyzed chemical oxidation at circumneutral pH under different conditions (H 2 O 2 or Na 2 S 2 O 8 with or without Fe(II) or magnetite as catalyst). Organic analyses were performed by using Gas Chromatography -flame ionization detector (GC-FID). The performance of chemical oxidation treatment in correlation with the availability of oil hydrocarbons in aged contaminated sediments was discussed.

. Chemical reagents and materials

  Following reagents were used: FeSO 4 .7H 2 O with purity greater than 98% obtained from VWR, sodium persulfate (Na 2 S 2 O 8 ) and hydrogen peroxide 35% (H 2 O 2 ) of Acros Organics.

Table 2

 2

	. These sediments were

Table 2 :

 2 Salient characteristics of tested sediments

	TMH	PIC

* represents the sum of clay, fine silt and coarse silt.

Acknowledgments:

Authors thank the GISFI (French Scientific Interest Group -Industrial Wasteland, www.gisfi.univ-lorraine.fr).