Serigne Gueye

Sophie Michel

Adnan Yassine

Solving a Berth Assignment Problem

Keywords: Vessels Assignment, Port Terminal, 0-1 Quadratic Programming, Linearization, Linear Programming, Constraint Programming One may consider

The Berth Allocation Problem (BAP) is the problem of allocating berthing spaces and scheduling container vessels on these spaces so as to minimize total weighted time. We study a version of BAP in which containers are moved between vessels and berth space is abundant. Thus, the problem reduces to optimally assign vessels to berths. We call it the Berth Assignment Problem (BASP). We formulate it as a non standard Quadratic Assignment Problem, and we show that BASP is NP-Hard. The formulation is simplified, linearized, and valid inequalities are found. Numerical results are shown.

Introduction

The Berth Allocation Problem (BAP) is the problem of allocating berth spaces to vessels in a container terminal. The general BAP problem stands as follows.

Along the wharf of a terminal, n container vessels, with known arrival dates a i , have to be berthed. Depending on their sizes, as well as the number of containers in the vessel, each of them needs a known length along the wharf. To operate (uploading and downloading operations) a vessel a processing time p i is necessary. The problem is then to find when (date u i ) and where a vessel has to be assigned in order to minimize some costs without exceeding the wharf capacity (length). Several Berth Allocation Problems have appeared in the litterature differing on the assumptions made.

anywhere, or as a collection of sections (i.e subdivisions) in which vessels are allowed to berth. The first case is known as the continuous case while the second is the discrete one. The costs to optimize is in general the vessel waiting times computed as follows min n i=1 w i (u i + p i -a i ) i.e sum of waiting times But other objectives may be considered. For instance, H. Wai [START_REF] Wai | Berth allocation and planning. (cover story)[END_REF] try to solve a multiobjective problem in which one of the objective is the waiting time and the other a container transhipment cost. In this type of problem, it is supposed that between vessels arriving in the port flows of containers exist. These flows are summarized in a flow matrix of containers giving how many containers a vessel has to transfer to another one and vice versa. For a complete survey of models and solution methods for the Berth Allocation Problem, the reader may be refered to Guan and Cheun [START_REF] Guan | The berth allocation problem : Models and solution methods[END_REF] , Cordeau et al [START_REF] Cordeau | Models and tabu search heuristics for the berth allocation problem[END_REF], Lim [START_REF] Lim | The berth planning problem[END_REF], Wang and Lim [START_REF] Wang | A stochastic beam search for the berth allocation problem[END_REF] among others.

In the problems described above, the operational details are omitted. Indeed, the way that containers are exactly handled are not taken into account. For instance, in transhipment of containers, a direct transhipment (from vessel to vessel) is not always necessary . A container can be uploaded on the yard area, then stocked and moved latter in another yard area in order to be downloaded in another vessel. These operational aspects are relaxed. As a consequence, allocating vessels on berths may be viewed as a tactical step of a port logistic management.

We are interested in another special case of BAP. This case was proposed by the Information Technology Company GTI (Le Havre) as a subproblem embedded in a Transport and Logistic package. We consider the discrete case in which the wharf is divided in sections of equal size. n vessels have to be moored and each of them needs a number b i (i = 1, 2, ..., n) of sections. We suppose that there are enough sections for all the vessels and that some containers have to be moved between vessels. As in the Wai problem the container flows are summarized in a matrix. Since the number of sections is sufficient, the time dimension dropped and the problem consists only in assigning optimally vessel to sections in such a way to minimize container transhipment costs. That is the reason why we call this problem a Berth Assignment Problem (BASP) instead of a Berth Allocation Problem (BAP). Because of the number of vessels to consider in practice (says ≥ 10) and the combinatorial complexity of the problem (n! possible assignments), solving such problems in a "reasonable" amount of time is only possible using heuristics. We have proposed in 2005 a greedy heuristic (see Gueye [START_REF] Gueye | Optimisation de l'affectation de "trades[END_REF]). Nevertheless, the evaluation of the solution quality is possible only in comparison with the optimal one's. Because lower is the gap between the heuristic solution value and the optimal value and better is the heuristic. In this paper, we are interested in solving exactly this problem.

The BASP is explained in detail in section 2. We modelize it in section 3 as a non standard version of the well known Quadratic Assignment Problem. This version is NP-Hard as the QAP. We propose a linearization of the model in the same section. This linearization leads theoretically to a very poor lower bound. To improve it some valid inequalities are proposed in section 4. The integer linear problem presents in section 3 and the valid inequalities are used together in an algorithm. This algorithm has been implemented and tested on instances derived from QAPLIB Benchmark. Numerical results are shown in section 5. Finally, some perspectives about this work are point out in section 6.

The Berth Assignment Problem

Let us define a section as an area, of known length, along the wharf, in which some containers may be downloaded or uploaded from a vessel. This definition is schematically illustrated in the figure below where the bold line represents a section. We define a "location" as a contiguous series of sections. The size of a location is the number of its sections. If n is the number of vessels (resp. locations) then each vessel (resp. location) may be identified by an index i = 1, 2, ..., n (resp. k = 1, 2, ..., n). If we consider a wharf as the horizontal line below, locations are arranged starting from the left-hand-side, to the right-hand-side, in increasing order of the indices. In this example the size of the location 1, 2, 3 are respectively 3, 2, 1.

Vessels arriving in the port have to be assigned in locations. For each vessel the length of its location depends on its size, but also on the number of containers to handle in the vessel. We will notice b i the location size of the vessel i. If S is the total number of sections of the wharf we must have

n i=1 b i = S.
It is important to observe that if k is the index of a location then its size depends on the vessel assigns to it. In other words, the size of a location is not known a priori but it is determined only when a vessel is assigned to the location. Therefore, let π be an assignment of vessels on locations, we will notice s(π, k) the size of the location k.

Vessels, assigned along the wharf, exchange some containers, mathematically expressed by a traffic matrix F = {f ij } i,j=1,2,...,n where i and j are vessel indices, and f ii = 0 for all i.

Let π be an assignment where π(i) = k (i.e vessel i is assigned to location k), π(j) = l (i.e vessel j is assigned to location l). The cost of moving f ij containers from k to l is f ij d kl (π), where d kl (π) is the distance (depending on the assignment π) between locations k and l. This distance is evaluated as follows.

Let us notice α (meters) the length of a section. For k < l, we have

d kl (π) = ( s(π, k) 2 + l-1 u=k+1 s(π, u) + s(π, l) 2 )α meters. (1)
The overall cost of the assignment π is then :

n i=1 n j=1 f ij d π(i)π(j) (π) meters.
Example 2.1 Let us consider the following example to illustrate mathematical expressions above.

In this example, 3 vessels (1, 2, 3) are assigned to 3 locations (1, 2, 3). This assignment corresponds to the identity permutation π(i) = i, i = 1, 2, 3.

The traffic matrix between vessels is given by

F =        0 10 20 8 0 2 10 10 0       
Let us consider vessels 1 and 3 assigned respectively to locations 1 and 3. Then, if the size of a section is 50 meters then the distance between locations 1 and 3 will be

d 13 (π) = 150 2 + 50 + 100 2 = 175 m
More generally, with the assignment π, we have the following distance matrix between locations : In this case, we have the following different distance matrix between locations : The movements of containers between locations need human and material resources, and are time-consuming. Higher is the distance between two vessels and higher will be the economic cost induced by the container movements. This imply that lower is the total overall distance to transfer all containers and lower is the associated economic cost for the port manager. Thus, we have to solve the following optimization problem :

D(π) =        0
D(π * ) =        0
(BASP ) : M in π∈Sn n i=1 n j=1 f ij d π(i)π(j) (π).
where S n is the permutation set of {1, 2, ..., n}, F = {f ij } i,j=1,2,...,n is a known traffic flow matrix and D(π) = {d kl (π)} k,l=1,2,...,n is a distance matrix whose values depend on the assignment π. This mathematical programming problem (BASP ) may be viewed as a "non standard" version of the well known Quadratic Assignment Problem (QAP). At the opposed of the "standard" QAP where the distance matrix D is fixed (i.e does not depends on any assignment π), in the BASP this matrix may change.

The BASP assignment is a very difficult combinatorial optimization problem (NP-Hard). Hence, because of the number of vessels to consider in practice and the complexity of the problem, we have proposed in 2005 (see [START_REF] Gueye | Optimisation de l'affectation de "trades[END_REF]) a greedy heuristic for the IT consulting company GTI. In the sequel, a formulation and an exact resolution scheme are now developped.

Mathematical formulation

Let us denote by i, j ∈ {1, 2, ..., n} the vessel indices and by k, l ∈ {1, 2, ..., n} the location indices. As in the standard QAP formulation, we define the following binary variables

x ik =      1 if vessel i is assigned to location k 0 otherwise
These variables have to verify the standard QAP constraints :

n i=1
x ik = 1 and n k=1

x ik = 1.

Let x = {x ik } 1≤i,k≤n be the variables matrix. Since an assignment of x corresponds to a permutation, the overall cost that has to be minimized in the problem is

n i=1 n j=1 n k=1 n l=1 f ij d kl (x)x ik x jl .
The expression of d kl (x) is obtained as follows. Let us first recall the equation ( 1) d kl (π) for k < l :

d kl (π) = ( s(π, k) 2 + l-1 u=k+1 s(π, u) + s(π, l) 2 )α meters,
where s(π, k) is the size of the location k in the case of permutation π. Since, in our notation, x stands now for π, we can also write bmx ml )α,

d kl (x) = ( s(x, k) 2 + l-1 u=k+1 s(x, u) + s(x, l) 2 )α meters.
1 ≤ k, l ≤ n, |k -l| ≥ 2 (2) d kl (x) = ( 1 2 n m=1 bmx mk + 1 2 n m=1 bmx ml )α, 1 ≤ k, l ≤ n, |k -l| = 1 (3) d kk (x) = 0, 1 ≤ k ≤ n (4) n i=1 x ik = 1, 1 ≤ k ≤ n (5) n k=1 x ik = 1, 1 ≤ i ≤ n x ∈ {0, 1} n 2
This is 0-1 quadratic problems for which we propose the following linearization.

Since α (the section length) is a constant, we can remove this value to constraints (1) and ( 2). Moreover, a formulation with only variables x, flow matrix F and sizes b can be derived by replacing d kl (x) with its corresponding explicit expressions (constraints (1),( 2),( 3)). Now, let y kml be the following variables :

y kml = max{k,l}-1 u=min{k,l}+1
x mu , k, l, m = 1, 2, ..., n and |k -l| ≥ 2 y kml x ik x jl . We have shown that using these variables BASP may be linearized as follows :

y kml = 0, k, l, m = 1,
(LBASP) Min n i=1 n j=1 i =j n m=1 m =i,j f ij bmt imj s-t : (4) n i=1 x ik = 1, 1 ≤ k ≤ n (5)
n k=1 x ik = 1, 1 ≤ i ≤ n (9) t imj ≥ max{k,l}-1 u=min{k,l}+1 xmu + x ik + x jl -2 1 ≤ i, j, m, k, l ≤ n, i = j, m = i, j, |k -l| ≥ 2 t imj ≥ 0 1 ≤ i, j, m ≤ n, i = j, m = i, j x ∈ {0, 1} n 2
This is the formulation that we seek to solve using linear integer programming techniques, in particular Branch-and-Bound algorithm. It is well-known that the corresponding linear relaxation gives a lower bound of the optimal value of the problem (BASP). It is also well-known that to reduce the Branch-and-Bound tree size, as well as the processing time, the gap (between the lower bound and the optimal value) should be as weak as possible. Unfortunately, the linear relaxation of LBASP gives the poorest bound that may be expected. Indeed, the optimal value of the relaxation is always equal to 0 as proved in the property below.

Propriety 3.1 Let (LBASP ) the linear relaxation of (LBASP ), and V (LBASP ) the optimal value of (LBASP ). We have V (LBASP ) = 0.

This result show that it remains to find many valid inequalities of the (LBASP ) feasible domain. The section 4 below deals with these inequalities.

Valid inequalities

We begin this section with the following property. It will be useful to reduce the number of constraints of (LBASP ), as well as to find some valid inequalities.

Propriety 4.1 .

(i) t imj = 1 if and only if vessel m is between vessel i and j in the optimal assignment.

(ii) t imj = t jmi Notice, for the point (i), that when t imj = 1 we only know that m is between i and j without any information about the position of i (resp. j) in comparison to m. In other words, we don't know if i (resp. j) is at the left (resp. right) of m or not. With this property it follows that only three cases are (exclusively) possibles : m is between i and j or i is between m and j, or j is between i and m. We thus derive our first valid equalities below Proposition 4.2 .

t imj + t mij + t ijm = 1 ∀ 1 ≤ i < m < j ≤ n.
Taking into account these equalities and the symmetry property (ii), (LBASP) may be strenghthened as follows

(LBASP * ) Min n i=1 2 n j=i+1 n m=1 m =i,j f ij bmt imj s-t : (4) n i=1 x ik = 1, 1 ≤ k ≤ n (5)
n k=1 x ik = 1, 1 ≤ i ≤ n (9) t imj ≥ max{k,l}-1 u=min{k,l}+1 xmu + x ik + x jl -2 1 ≤ i, j, m, k, l ≤ n, i < j, m = i, j, |k -l| ≥ 2 (10) t imj ≥ 0 1 ≤ i, j, m ≤ n, i = j, m = i, j (11)
t imj + t mij + t ijm = 1 ∀ 1 ≤ i < m < j ≤ n x ∈ {0, 1} n 2
Let observe that the linear relaxation of this formulation discard the solution x ik = 1 n , t imj = 0. However, constraints (11) do not described the entire convex hull of (LBASP ). We have found an exponential number of other valid inequalities, derived from the meaning of variables t imj (see property 4.1 (i)), and by observing the optimal solution of the linear relaxation. Because of the huge number of valid inequalities, it is not possible, for computer memory reasons, to add all of them simultaneously. Inequalities have to be added progressively, using an iterative algorithm for which at each iteration a so-called "separation problem" is solved and gives as a result a new inequality for the formulation. Such type of algorithm is called a cut (or constraint) generation algorithm in the litterature. We have implemented this type of algorithm. At each iteration the separation problem is solved using Constraint Programming (CP) Techniques. To obtain accurate details about the algorithm, the separation problem and CP techniques used we refer the reader to the technical report (Gueye [4]).

Numerical Tests

With our algorithms, some numerical tests have been performed. Since no practical datas or academic benchmarks were available , we generate ourselves some instances. As the Berth Assignment Problem is a variant of the Quadratic Assignment Problem, our instances have been derived from QAPLIB benchmark. Each instance is define by the container flow matrix, noticed F , and by the number of sections occupied by each vessel, noticed b i (i = 1, 2, ..., n).

Numerical experiments have been performed on a Personal Computer Dell, with 2 G0 RAM, Intel Core Processor of 1.86 GHz on linux Suze Operating System. Ilog Cplex 9.3 and Solver 9.3 have been used. The results are reported in table 1. We indicate in column 3 of this table the optimal value of (LBASP * ) (i.e the linear relaxation of (LBASP * )), and in colum 4 the optimal value of (BASP ). In column 5 the number of nodes of the Cplex 9.3 branch and bound tree is reported. We have fixed the Cplex parameters MIPEmphasis to 1 in order to work toward a rapid feasibility. The experiences show that, in our formulation, emphasizing on feasibility is better (less computational time) than emphasizing on optimality. The Cplex MIP heuristic frequency has been switched off. The processing time has been limited to 1 hour. When the processing time exceed this value, no result is given in the table. All other parameters have been leaved at their default values.

The maximal number of vessels that we are able to assign optimally is currently lower or equal to 10. In all cases, the gap between lower bounds and optimal values do not exceed 12% from the optimal values. But, the computational times grow exponentially with the size n. This is explain by the huge amount of constraints, in particular constraints 9, of the formulation (LBASP * ). The instances chr12a.dat and chr12b.dat seem easiers since no gap exists.

Conclusion

This paper deals with a difficult combinatorial optimization problem. In this problem vessels have to be assigned along the wharf of a port terminal in such a way to minimize the total transhipment cost of containers. The problem has been formulated as a non standard version of the Quadratic Assignment Problem. It has been shown that this version is NP-Hard. Following many reformulation steps, and a linearization, we have obtained a integer linear formulation. But this formulation leads to poor lower bound. Without adding new valid inequalities, it is not possible to solve optimally the problem within a reasonable amount of time. We have found an exponential number of such valid inequalities. Then, a cut generation algorithm has been implemented to generate these inequalities. At each iteration of our alogrithm, Constraint Programming techniques are used to find the cuts.

Currently, the maximal number of vessels that we may assign optimally is 10. This means that it remains further works in order to improve the numerical results. Such improvement may be reach considering another formulation with less constraints than the one presented in this paper.

Following the improvements, we seek also to take into account the arrival time of the vessels. This aspect will lead to a challenging scheduling problem and will make the problem formulation closer to the reality.

  And the overall cost of π is (10+8) * 100 + (20+10) * 175 + (10+2) * 75 = 7950 meters (i.e 7 km 950 m) Let us observe what happens with the other assignment π * (1) = 1, π * (2) = 3, π * (3) = 2 shown below.

  And the overall cost of π * is (10+8) * 200 + (20+10) * 125 + (10+2) * 75 = 8250 meters (i.e 8 km 250 m) 2

bbf

  By definition of s(x, k), x and b i (number of sections need for vessel i), we haves(x, k) = n m=1 b m x mk . Thus, for k < l, d kl (x) is d kl (x) = ( m x ml )α.Hence, for any values for k and l with k = l d kl (x) m x ml )α. When k = l, we have d kl (x) = 0. The notations above give the complete mathematical program of the Berth Assignment Problem (BASP) : ij d kl (x)x ik x jl s-t : (1) d kl (x) = (

2 ,

 2 ..., n and |k -l| ≤ 1, and t imj =


Aliou Baldé
email: aliou.balde@iplesp.upmc.fr

Laurence Lièvre

Almoustapha Issiaka

Fodié Diallo

Issouf Alassane Maiga

Dominique Costagliola

Sophie Abgrall

Risk factors for loss to follow-up, transfer or death among people living with HIV on first ART regimen in Mali

Keywords: HIV, access to care, risk factors, loss to follow-up, Mali

Objectives: Risk factors for loss-to-follow-up (LTFU) were assessed for people living with HIV (PLHIV) at various reference outpatient clinics (expertise level II) and hospitals (expertise level III) in Mali.

Methods

: HIV-1-positive adults starting antiretroviral therapy (ART) in 2006-2013 were eligible. Risk factors for LTFU, defined as no visit in the six months preceding the last database update, were assessed with Cox model, taking into account the competing risks of transfer and death. Potential risk factors at start of ART were demographic and socioeconomic variables, WHO stage and CD4, period, ART regimen, region of care, expertise level and distance from home.

Results: We included 9,821 PLHIV, 33% male, starting ART at nine outpatient clinics and seven hospitals (5+2 in Bamako and 4+5 in the regions) with a median (IQR) CD4 count of 153/µL (56-270). Five-year cumulative incidences of LTFU, transfer and death were 35.2%, 9.7% and 6.7%, respectively. People followed at Bamako hospitals > 5 km from home, at regional hospitals or at regional outpatient clinics < 5 km from home were at higher risk of LTFU than people followed at Bamako outpatient clinics, whereas people followed at regional outpatient clinics 5-50 km away from home were at lower risk for LTFU.

Deaths were less frequent at hospitals, whether in Bamako or the regions, than at Bamako outpatient clinics, and more frequent at regional outpatient clinics.

Conclusions:

Expertise level and distance to care were associated with LTFU. Stigmatisation may play a role for PLHIV living close to the centres in the regions.

Introduction

Loss to follow-up (LTFU) is a predictor of poor survival in PLHIV, and is the principal problem in HIV care and in the evaluation of ART access programmes in sub-Saharan Africa [START_REF] Ekouevi | International epidemiological Database to Evaluate AIDS (IeDEA) West Africa Collaboration. Low retention of HIV-infected patients on antiretroviral therapy in 11 clinical centres in West Africa[END_REF][START_REF] Cornell | International epidemiological Database to Evaluate AIDS (IeDEA) Southern Africa Collaboration. Temporal changes in programme outcomes among adult patients initiating antiretroviral therapy across South Africa, 2002-2007[END_REF][START_REF] Grimsrud | Outcomes of antiretroviral therapy over a 10-year period of expansion: a multicohort analysis of African and Asian HIV programs[END_REF], in which mean attrition rates of 25.0% and 38.4% at 12 and 24 months have been reported respectively, with LTFU accounting for 56% of attrition [START_REF] Rosen | Patient retention in antiretroviral therapy programs in sub-Saharan Africa: a systematic review[END_REF].

In Mali, a West African country with 16.8 million inhabitants, a non-governmental organisation, "Association de Recherche, de Communication et d'Accompagnement à Domicile des personnes vivant avec le VIH/Sida" (ARCAD-SIDA), began to provide healthcare for PLHIV in Bamako (the capital of Mali) in 1998, gradually expanding its activities to several regions after November 2001. Healthcare for PLHIV became free after mid-2004, with support from the Global Fund, "Ensemble pour une Solidarité Thérapeutique Hospitalière en Réseau" (ESTHER), a French agency for AIDS care in developing countries, and the Malian government [START_REF]Cellule Sectorielle de Lutte contre le VIH-Sida du Ministère de la Santé et de l'Hygiène Publique[END_REF][START_REF] Raguin | The ESTHER hospital partnership initiative: a powerful levy for building capacities to combat the HIV pandemic in low-resource countries[END_REF]. There were 99,000 adults living with HIV in Mali in 2013, 27,000 of whom were treated [START_REF]Cellule Sectorielle de Lutte contre le VIH-Sida du Ministère de la Santé et de l'Hygiène Publique[END_REF].

Our aim was to analyse ART outcomes in ESTHER-supported clinical centres after more than ten years of the scaling-up of access to ART, and to assess risk factors for LTFU, transfer to other centres and death, focusing on structural barriers, such as care in the capital or in the regions, care centre expertise level, distance from home to the centre.

Methods

Study setting

In Mali, HIV care and ART are provided by community outpatient clinics (n=14) offering consultations with nurses and/or general practitioners (expertise level I), reference outpatient clinics (n=54) offering consultations with general practitioners, obstetricians and ophthalmologists, minimal medical tests, including CD4, and simple radiology (level II), or hospitals (n=9) offering the same care plus specialized consultations and medical tests, including viral load (VL), and CT scans (level III).

The ESTHER initiative focused on reinforcing the clinical skills of health professionals, simultaneously providing clinical mentoring, teaching, monitoring of CD4, VL and HIV drug resistance. In addition, since 2005, ESTHER has provided an electronic medical records system "Evaluation et Suivi Opérationnel des Programmes d'ESTHER" (ESOPE, Epiconcept, France) for monitoring PLHIV care. This system was progressively implemented at five reference outpatient clinics and two hospitals in Bamako, and at four reference outpatient clinics and five hospitals in five regions (Supplementary figure 1) [START_REF] Raguin | The ESTHER hospital partnership initiative: a powerful levy for building capacities to combat the HIV pandemic in low-resource countries[END_REF]. The PLHIV included in the ESOPE database accounted for 64% of the 51,000 PLHIV in healthcare for HIV in Mali at the last database update, on March 31 st 2015 [START_REF]Cellule Sectorielle de Lutte contre le VIH-Sida du Ministère de la Santé et de l'Hygiène Publique[END_REF].

Study population

This study included all the HIV-1-positive adults, aged ≥18 years who began ART between January 1 st 2006 and December 31 st 2013 at one of the 16 care centres contributing to the ESOPE database, thus beginning treatment at least one year before the last database update (31 March 2015), and who returned for their one-month visit. We did not include PLHIV beginning ART before 2006, because the first national guidelines on HIV care and ART were provided at the end of 2005, with subsequent revisions (2008, 2010 and 2013) according to WHO guidelines. Treatment was initially supplied for one month, then provided for three months at a time, with clinical visits scheduled every three months and biological tests performed every six months [START_REF]Cellule Sectorielle de Lutte contre le VIH-Sida du Ministère de la Santé et de l'Hygiène Publique[END_REF]. In 2005, guidelines recommended treatment for all PLHIV with CD4 ≤ 200 cells/µL, or at WHO stage 3 with CD4 ≤ 350 cells/µL or at WHO stage 4. In 2013, ART was recommended for all PLHIV with CD4 ≤ 500 cells/µL or at WHO stage 3-4 [START_REF]Cellule Sectorielle de Lutte contre le VIH-Sida du Ministère de la Santé et de l'Hygiène Publique[END_REF].

Definitions

The outcome measures were LTFU, transfer to another care centre or death. LTFU was defined as an interval of at least six months between the last clinical visit and the last database update at each care centre, for PLHIV not known to be dead or to have been transferred elsewhere, because it is the definition yielding the lowest rate of misclassification [START_REF] Chi | International Epidemiologic Databases to Evaluate AIDS (IeDEA) Collaboration[END_REF]. As clinicians first consider an individual to be LTFU at the time of the first missed visit, date of LTFU was the date of the last centre visit plus three months [9-10].

Statistical analyses

The follow-up of individuals began with the initiation of ART and was censored at the first of the three competing events (LTFU, transfer or death), after five years on ART or in March 2015, whichever occurred first. We estimated the cumulative incidence of these three events over time and associations of individual characteristics at the start of ART with these events, using univariable and multivariable competing risk regression models [11].

The characteristics of PLHIV recorded at the start of ART were: sex, presence or absence of pregnancy, age, WHO stage and CD4, period of ART initiation, ART regimen, a variable combining the location (Bamako or the regions), expertise level (reference outpatient clinic or hospital) of the care centre and its distance from the individual's home, marital status, level of education and professional activity. The ART initiation periods were selected to coincide with the changes in national guidelines over time. The distance from home to the care centre was determined, in kilometres (km), with Google Maps, by calculating the distance, by road, between the individual's village or quarter of residence and his care centre.

Multiple imputations were performed with chained equations [12] to deal with missing values, assuming data were conditionally missing at random. More specifically, we fitted linear regression, logistic regression and multinomial regression models for CD4 (fourth root), WHO stage, and educational level, respectively. All models included age, sex, pregnancy, period of ART initiation, region of care and care expertise level, marital status, death and the follow-up time in days (log 10 ) as covariables. We created ten imputed datasets; analyses were run separately on each dataset and the results were combined by Rubin's method [13].

Analyses were conducted with SAS 9.4 (SAS Institute Inc., NC, USA) and STATA 12.1 (STATA Corporation, College Station, TX, USA). The study was approved by the national AIDS programme in Mali in collaboration with the non-governmental organisations ARCAD-SIDA and ESTHER. Written or oral informed consent was obtained from all PLHIV.

Results

Of the 17,730 HIV-positive adults contributing to the ESOPE database and starting ART between 2006 and 2013, 9,821 individuals with CD4 or WHO stage data available at start of ART were included in this study (Supplementary figure 2). Table 1 shows the characteristics of the PLHIV, overall and according to the region of care and the care centre expertise level at the start of ART. Median follow-up time was 29 months [interquartile range (IQR), 10-59]. Individuals with the worst prognostic factors (WHO stage 3-4 and CD4 < 200 cells/µL) were more frequent at regional outpatient clinics (50.7%) and less frequent at Bamako outpatient clinics (22.4%) (Table 1).

Incidence of LTFU, transfer and death

The cumulative incidence of LTFU rose from 15.7% (95% CI: 15.1-16.4) at one year to 35.2% (95% CI: 34.2-36.3) at five years on ART, of transfer to another centre rose from 5.1% (4.7-5.6) at one year to 9.7% (9.1-10.4) at five years, and of mortality rose from 4.0% (3.7-4.4) at one year to 6.7% (6.2-7.2) at five years (Supplementary figure 3).

Risk factors for LTFU and death

Five years after ART initiation, all of the variables studied were associated with LTFU in both univariable and multivariable analyses (Supplementary table 1). Compared to PLHIV followed at Bamako outpatient clinics, risk of being LTFU was approximately 45% higher and three times higher in PLHIV followed at Bamako hospitals if the hospital was located more than 5 km from the individual's home or at the regions hospitals regardless of distance respectively, and 29% higher and 19% lower in people followed at regional outpatient clinics and living less than 5 km from their care centre or living 5-50 km respectively (Figure 1A). Deaths were less frequent at hospitals, whether in Bamako or the regions, than at Bamako outpatient clinics, and more frequent at regional outpatient clinics (Figure 1A).

Figure 1B shows the adjusted risk factors for LTFU by distance from home for each type of care centre. LTFU rates did not differ according to the distance between home and the care centre among the PLHIV followed at Bamako outpatient clinics. People treated at Bamako hospitals and living 5-50 km or more than 50 km away from their care centre were more likely to be LTFU than those living less than 5 km. For PLHIV treated at regional outpatient clinics, those living 5-50 km or more than 50 km away from their care centres were less likely to be LTFU than those living less than 5 km (Figure 1B).

Discussion

Between 2006 and 2015 in Mali, 15.7% and 35.2% of PLHIV were LTFU within one and five years of starting ART respectively. People followed at Bamako hospitals 5 km or more from their homes, at regional hospitals or at regional outpatient clinics less than 5 km from their homes were at higher risk of LTFU than those followed at Bamako outpatient clinics, whereas people followed at regional outpatient clinics 5-50 km from their homes were at lower risk of LTFU. Estimated death rates were low, whereas LTFU rates were relatively high. There is no registry of deaths in Mali and the national programme has no dedicated resources for the active tracing of PLHIV [START_REF] Cornell | International epidemiological Database to Evaluate AIDS (IeDEA) Southern Africa Collaboration. Temporal changes in programme outcomes among adult patients initiating antiretroviral therapy across South Africa, 2002-2007[END_REF]14].

We were therefore unable to estimate the true mortality rate. It is also possible that mortality rates are lower for the more recent periods in which universal ART has been recommended [15]. Moreover, the main objective of the study was to estimate LTFU before death whilst taking deaths among PLHIV remaining in healthcare system into account. Transfers to another care centre were also taken into account as competing risks, to improve the estimation of LTFU and death rates. The main strength of this study is the large number of PLHIV included, corresponding to 64% of those in healthcare in Mali and followed in either the capital or the regions. Without active tracing of people LTFU, Cornell et al [START_REF] Cornell | International epidemiological Database to Evaluate AIDS (IeDEA) Southern Africa Collaboration. Temporal changes in programme outcomes among adult patients initiating antiretroviral therapy across South Africa, 2002-2007[END_REF] observed in South Africa a risk of LTFU of 14.4% (95% CI, 14.1-14.8) and 35.8% (95% CI, 34.3-37.4) at 1 and 5 years respectively, those estimates being close to ours.

Distance from home to care centres is a major obstacle to access to care, and many sub-Saharan countries have tried to decentralise HIV care from hospitals to healthcare centres and other primary health facilities closer to the community [16]. Decentralisation is associated with a shorter travelling distance, higher rates of enrolment in ART programmes and higher levels of retention in healthcare, even if the decrease in distance is only moderate [17]. We also showed that living more than 5 km from the healthcare centre was associated with a higher risk of LTFU for the hospitals, whether in Bamako or the regions, than for Bamako outpatient clinics. Objective measures of transportation are, therefore, required, to optimise PLHIV access to HIV care in resource-limited settings [START_REF] Siedner | GPS-measured distance to clinic, but not self-reported transportation factors, are associated with missed HIV clinic visits in rural Uganda[END_REF].

At outpatient clinics, no difference in LTFU with distance was observed in Bamako, whereas, in the regions, people living less than 5 km from the healthcare centre were more likely to be LTFU than those living more than 5 km away, raising questions about the possibility of HIV-related stigma and secrecy [START_REF] Hurley | The role of patient-provider communication in engagement and re-engagement in HIV treatment in Bamako, Mali: a qualitative study[END_REF]. Better patient-provider communication is associated with a better engagement and re-engagement in HIV, as disrespect, discrimination and scolding are the major reasons for PLHIV disengaging from ART [START_REF] Hurley | The role of patient-provider communication in engagement and re-engagement in HIV treatment in Bamako, Mali: a qualitative study[END_REF].

It is noteworthy that distance was not associated with death in our study as previously reported [START_REF] Bassett | Barriers to care and 1year mortality among newly diagnosed HIV-infected people in Durban, South Africa[END_REF], but death rates were lower in the hospitals than at outpatient clinics, whatever the distance, highlighting the need for improvements in healthcare centre expertise.

We identified several structural risk factors associated with LTFU or death, some of which are modifiable, such as the number and expertise of care centres. Possible stigmatisation for individuals travelling very short distances to the care centres in the regions and insufficient care provision for individuals living at much greater distances should be investigated. ART programmes may increase the success with which PLHIV are retained in the healthcare system, by investing in skills and resources for individual problemsolving rather than offering the standard "education" [START_REF] Hurley | The role of patient-provider communication in engagement and re-engagement in HIV treatment in Bamako, Mali: a qualitative study[END_REF]. A knowledge of the risk factors associated with programme attrition, LTFU or death would help clinicians to identify PLHIV requiring additional support to remain in the healthcare system, and the opening of other healthcare facilities would contribute to the achievement of the UNAIDS targets 90-90-90. Data are counts (proportions) and medians (interquartile range). The percentages of imputed missing data were 15% for CD4 cell count, 8% for WHO stage and 16% for educational level Abbreviations: SHR, subdistribution hazards ratio; CI, confidence interval; NNRTI, non-nucleoside reverse transcriptase inhibitor; PI, protease inhibitor; EFV, efavirenz; NVP, nevirapine; NRTI, nucleoside reverse transcriptase inhibitor, AZT, zidovudine; d4T, stavudine; 3TC, lamivudine; TDF, tenofovir; FTC, emtricitabine. a Others: PI (lopinavir/ritonavir, indinavir/ritonavir, saquinavir/ritonavir and atazanavir/ritonavir) or NRTIs in 3-NRTI combinations. b Others NRTI: abacavir (ABC)/didanosine (ddI), ABC/3TC, ddI/3TC, TDF/ddI and TDF/ABC. c Married or living with a partner. d Other for marital status: single, divorced, widowed and missing (6%). e Professional activities: professional activity categories combined types and sectors of employment, housewife without any professional activity, public or private sector salaried employment, various non-salaried activities, such as being self-employed, farmer/fisherman and other (unemployed and missing (3%)). *SHR adjusted for variables at the start of ART: sex and pregnancy, age, WHO stage and CD4 count, period of ART initiation, initial ART regimen, marital status, educational level and professional activities. Abbreviations: SHR, subdistribution hazards ratio; CI, confidence interval; NNRTI, non-nucleoside reverse transcriptase inhibitor; PI, protease inhibitor; EFV, efavirenz; NVP, nevirapine; NRTI, nucleoside reverse transcriptase inhibitor, AZT, zidovudine; d4T, stavudine; 3TC, lamivudine; TDF, tenofovir; FTC, emtricitabine. e Professional activities: professional activity categories combined types and sectors of employment, housewife without any professional activity, public or private sector salaried employment, various non-salaried activities, such as being self-employed, farmer/fisherman and other (unemployed and missing (3%)).

follow-up in HIV treatment programs: a statistical analysis of 111 facilities in

  Africa, Asia, and Latin America. PLoS Med 2011; 8:e1001111. 9. Mberi MN, Kuonza LR, Dube NM, Nattey C, Manda S, Summers R. Determinants of loss to follow-up in patients on antiretroviral treatment, South Africa, 2004-2012: a cohort study. BMC Health Serv Res 2015; 15:259. 10. Shepherd BE, Blevins M, Vaz LM, Moon TD, Kipp AM, José E et al. Impact of definitions of loss to follow-up on estimates of retention, disease progression, and mortality: application to an HIV program in Mozambique. Am J Epidemiol 2013; 178:819-828. 11. Putter H, Fiocco M, Geskus RB. Tutorial in biostatistics: competing risks and multistate models. Statist Med 2007; 26:2389-2430. 12. White IR, Royston P, Wood AM. Multiple imputation using chained equations: Issues and guidance for practice. Stat Med 2011; 30: 377-399. 13. Rubin DB. Multiple imputation for nonresponse in surveys. New York: Wiley, 1987. 14. Egger M, Spycher BD, Sidle J, Weigel R, Geng EH, Fox MP et al; International epidemiological Database to Evaluate AIDS (IeDEA) East Africa, West Africa and Southern Africa. Correcting mortality for loss to follow-up: a nomogram applied to antiretroviral treatment programmes in sub-Saharan Africa. PLoS Med 2011; 8:e1000390. 15. Geng EH, Odeny TA, Lyamuya RE, Nakiwogga-Muwanga A, Diero L, Bwana M et al. Estimation of mortality among HIV-infected people on antiretroviral treatment in East Africa: a samplingbased approach in an observational, multisite, cohort study. Lancet HIV 2015; 2:e107-116. 16. Kredo T, Ford N, Adeniyi FB, Garner P. Decentralising HIV treatment in lower-and middleincome countries. Cochrane Database Syst Rev 2013; 6:CD009987.

  Abbreviations. LTFU, loss to follow-up; SHR, subdistribution hazards ratio; CI, confidence interval.

Figure 1 .

 1 Figure 1. A) Adjusted subdistribution hazard ratios (SHRs) for loss to follow-up (LTFU), transfer to another centre and death according to a variable combining regions of care, care expertise level and distance (km) from home and B) Adjusted SHRs for LTFU taking into account transfer to another centre and death as competing risk events, for each type of care centre according to distance (km) from home (Cox model).

  Abbreviations: SHR, subdistribution hazards ratio; CI, confidence interval; NNRTI, non-nucleoside reverse transcriptase inhibitor; PI, protease inhibitor; EFV, efavirenz; NVP, nevirapine; NRTI, nucleoside reverse transcriptase inhibitor, AZT, zidovudine; d4T, stavudine; 3TC, lamivudine; TDF, tenofovir; FTC, emtricitabine. a Other: PI (lopinavir/ritonavir, indinavir/ritonavir, saquinavir/ritonavir and atazanavir/ritonavir) or NRTIs in 3-NRTI combination. b Other NRTI: abacavir (ABC)/didanosine (ddI), ABC/3TC, ddI/3TC, TDF/ddI and TDF/ABC. c Married or living with a partner. d Other marital status: single, divorced, widowed and missing (6%).

Table 1 .

 1 Characteristics of individuals, overall and by region and care centre type at the start of antiretroviral therapy (ART)

	Characteristics at start of ART	Overall	Bamako	Bamako	Regional	Regional
						outpatient clinics	hospitals	outpatient clinics	hospitals
					(n=9 821)	(n=5 069)	(n=2 091)	(n=1 367)	(n=1 294)
	Sex and pregnancy
	Men				3 281 (33.4)	1 648 (32.5)	747 (35.7)	426 (31.2)	460 (35.5)
	Non-pregnant women	5 734 (58.4)	2 902 (57.3)	1 192 (57.0)	848 (62.0)	792 (61.2)
	Pregnant women		806 (8.2)	519 (10.2)	152 (7.3)	93 (6.8)	42 (3.3)
	Age (years)			36 (29-43)	35 (29-43)	37 (30-45)	35 (29-43)	36 (30-43)
	<30				2 623 (26.7)	1 393 (27.5)	523 (25.0)	371 (27.1)	336 (26.0)
	30-40				3 709 (37.8)	1 946 (38.4)	768 (36.7)	505 (37.0)	490 (37.9)
	≥40				3 489 (35.5)	1 730 (34.1)	800 (38.3)	491 (35.9)	468 (36.1)
	CD4 count cell/µL		153 (56-270)	162 (58-274)	152 (59-271)	146 (55-270)	129 (43-254)
	WHO stage and
	CD4 count cell/µL
	3-4 and CD4<200		3 115 (31.7)	1 137 (22.4)	792 (37.8)	693 (50.7)	493 (38.1)
	3-4 and CD4≥200		1 482 (15.1)	532 (10.5)	384 (18.4)	362 (26.5)	204 (15.7)
	1-2 and CD4<200		2 830 (28.8)	1 841 (36.3)	480 (23.0)	146 (10.7)	363 (28.0)
	1-2 and CD4≥200		2 394 (24.4)	1 559 (30.8)	435 (20.8)	166 (12.1)	234 (18.1)
	Periods of ART initiation
	2006-2007			1 997 (20.3)	1 238 (24.4)	354 (16.9)	244 (17.9)	161 (12.5)
	2008-2009			1 877 (19.1)	820 (16.2)	480 (23.0)	279 (20.4)	298 (23.0)
	2010-2012			4 163 (42.4)	2 056 (40.6)	871 (41.6)	602 (44.0)	634 (49.0)
	2013				1 784 (18.2)	955 (18.8)	386 (18.5)	242 (17.7)	201 (15.5)
	First NNRTI/PI
	NVP				6 644 (67.7)	3 408 (67.2)	1 276 (61.0)	1 076 (78.7)	884 (68.3)
	EFV Other	a			2 782 (28.3) 395 (4.0)	1 428 (28.2) 233 (4.6)	730 (34.9) 85 (4.1)	262 (19.2) 29 (2.1)	362 (28.0) 48 (3.7)
	First NRTI backbone
	AZT or d4T plus 3TC		7 214 (73.5)	3 499 (69.0)	1 524 (72.9)	1 135 (83.0)	1 056 (81.6)
	TDF plus FTC or 3TC Other b		2 330 (23.7) 277 (2.8)	1 383 (27.3) 187 (3.7)	498 (23.8) 69 (3.3)	226 (16.5) 6 (0.5)	223 (17.2) 15 (1.2)
	Duration on ART (months)	29 (10-59)	32 (12-64)	31 (14-62)	24 (7-51)	20 (5-45)
	Marital status Monogamous	c		4 223 (43.0)	2 162 (42.6)	814 (38.9)	554 (40.5)	693 (53.6)
	Polygamous Other d			2 109 (21.5) 3 489 (35.5)	1 216 (24.0) 1 691 (33.4)	432 (20.7) 845 (40.4)	317 (23.2) 496 (36.3)	144 (11.1) 457 (35.3)
	Education level
	None				4 543 (46.2)	2 414 (47.6)	659 (31.5)	736 (53.8)	734 (56.7)
	Primary/Koranic		2 982 (30.4)	1 493 (29.5)	804 (38.4)	400 (29.3)	285 (22.0)
	Secondary			1 707 (17.4)	904 (17.8)	409 (19.6)	195 (14.3)	199 (15.4)
	Higher			589 (6.0)	258 (5.1)	219 (10.5)	36 (2.6)	76 (5.9)
	Professional activities	e
	Housewife			3 984 (40.6)	2 154 (42.5)	735 (35.1)	579 (42.3)	516 (39.9)
	Public sector staff		594 (6.0)	268 (5.3)	192 (9.2)	46 (3.4)	88 (6.8)
	Private sector staff		1 866 (19.0)	1 109 (21.9)	577 (27.6)	90 (6.6)	90 (6.9)
	Self-employed		1 588 (16.2)	764 (15.1)	282 (13.5)	349 (25.5)	193 (14.9)
	Farmer/Fisherman		605 (6.2)	255 (5.0)	89 (4.3)	124 (9.1)	137 (10.6)
	Other				1 184 (12.0)	519 (10.2)	216 (10.3)	179 (13.1)	270 (20.9)
	Distance to care centre (km)
	<5				2 412 (24.6)	1 404 (27.7)	350 (16.7)	654 (47.8)	4 (0.3)
	5-50				6 269 (63.8)	3 406 (67.2)	1 472 (70.4)	432 (31.6)	959 (74.1)
	≥50				1 140 (11.6)	259 (5.1)	269 (12.9)	281 (20.6)	331 (25.6)

Acknowledgements

We thank all the participants from Malian ART programme and research assistants from Mali and France. Conceived and designed the study: AB AIM DC SA. Analyzed the data: AB LL DC SA. Data provision: AB AIM FD IAM DC SA. Had full access to all the data and takes the responsibility for the integrity of the data and the accuracy of the data analysis: AB DC SA. Statistical expertise: DC SA. Drafted the article: AB DC SA. Interpretation of the data, critical revision of the article for important intellectual content and final approval of the article: AB LL AIM FD IAM DC SA. Aliou BALDE is supported by a doctoral contract from Sorbonne Université (ex Université Pierre et Marie-Curie, UPMC, Paris 6).

This research did not receive any specific grant from funding agencies in the public, commercial, or notfor-profit sectors.


Mehdi Daanoune

Raphaël Clerc
email: raphael.clerc@institutoptique.fr.

Bruno Flament

Lionel Hirsch

Physics of Trap Assisted Photo-Multiplication in Vertical Organic Photoresistors

Several experimental groups have reported recently an intriguing high level of gain (Photo-Multiplication) in vertical organic photo-resistance (as well as in other technologies, such as perovskite for instance). This mechanism is sometimes named as "Trap-Assisted Photo-Multiplication". This paper investigates the origin of this mechanism by the mean of drift diffusion simulations, analytical theory and experiments, considering the particular case of PCDTBT:PC 60 BM photoresistors, although some conclusions are likely to apply in other technologies.

It turns out that an excess of charges (induced by electron-holes carrier generation) may trigger additional carrier injection, leading to photomultiplication, under specific circumstances. We call this mechanism "gain by injection enhancement". Electron (resp. hole) trapping for P only (resp. N only) device can play this role efficiently. As these additional carriers came from contacts, significant dark current injection is thus needed to achieve large value of gain, explaining why this mechanism can occur only in P (or N) only photoresistors (and not photodiodes or intrinsic photoresistor, i.e. with mid gap contacts). In such devices however, the detectivity remains intrinsically limited by the high level of dark injection currents required to get gain, and consequently this type of device may be interesting in particular technologies where it is not possible to achieve low dark currents using photodiodes. However, penalized by the slow traps dynamic, the cut off frequency of these devices remains extremely low (< 100 Hz). Also, this gain takes high value only at low irradiance, making photoresistor responsivity light dependent.

All these results bring new light in analysing and optimising photoresistor, opening a large field of investigation to take advantage of gain by injection enhancement.

1) Introduction

Organic semiconductors are attractive materials for optoelectronics applications. One of the major asset of this technology is that they can be printed onto flexible, large area substrates at room temperature -a process that is both low cost and highly scalable for organic solar cells 1,2 or photodiodes arrays 3,4,5 . In particular, organic photodiodes based on the bulk heterojunction concept 6 , thanks to the use of printing processes, allows the realization of large area imagers 7 , potentially flexible, at reasonable cost, featuring detectivity comparable to silicon detectors in the visible range 8 . In the meantime, the reliability and the cut-off frequency, the typical weak points of organic photodiodes, has also drastically been improved 8 . All these progresses have made possible to design and produce innovative devices for specific applications in the field of X ray imaging 4 , human-machine interaction, and fingerprint recognition for instance. This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

To date, most of organic photodetectors are essentially vertical photodiodes architecture where the active layer is sandwiched between asymmetric electrodes with different work functions 9 . I-V characteristics are then obtained without intentional doping like in classical PN or PIN silicon photodiodes. One can also use two electrodes with similar work functions to perform organic photoconductors 3 (also called photoresistors or photo conductive cells). These devices have the combined advantages to potentially feature gain (also called photo-multiplication), i.e. to collect more than one electron for one photon absorbed. Nevertheless, this architecture give rise to larger dark current and lower cut-off frequency 3 . Cadmium-Sulfide photoconductors are for instance currently extensively used for low cost IR detection [START_REF] Vollmer | Infrared Thermal Imaging: Fundamentals, Research and Applications[END_REF] .

Recently, a significant number of experimental papers [START_REF] Nie | [END_REF]12 have reported photoconductors featuring a huge value of the gain (record value ≈ 1000 [START_REF] Nie | [END_REF] ) in organic based devices as well in as other materials, such as perovskite 13,14,15 for instance. These results are intriguing for several reasons. Firstly, in most papers, a large gain is obtained typically when defects, potentially operating like traps, are intentionally introduced. This mechanism is referred sometimes as "Trap-Assisted Photo-multiplication" 16 17 18 . These defects can be either C 60 19 , nano-particles such as ZnO 20 , PbSe 19 , PbS 21, 22 23 , CdSe 24 , CdTe 25 , variation of the polymer/PCBM blend ratio 18,26,27 , or interface charge trapping layer 12 . This observation may seem counterintuitive, as trapping or recombination, penalizing the collection of photo-generated carriers, are usually avoided in photo-detection. Secondly, as simulations have not been used so far to confirm the value of the measured gain, its physical origin remains unclear and subject to speculation. For instance, some authors attribute the gain to the presence of volume traps 19 , other to interface quality 28 , combining either interface traps or injection controlled by tunneling 12 . Moreover, the gain has been observed essentially in vertical structures (as illustrated in Figure 1), while conventional crystalline photoconductors are planar, and collect the charges by lateral contacts [START_REF] Sze | Physics of Semiconductor Devices[END_REF] . In consequence, the classical theory of photo-conductance [START_REF] Sze | Physics of Semiconductor Devices[END_REF][START_REF] Bube | Photoconductivity of Solids[END_REF] (which treats separately the transport of photo-generated carriers in the horizontal direction and the carrier generation by light in the vertical direction) is not supposed to apply in vertical device. In fact, even if this theory is sometimes claimed to explain experiments [19][20][21][22][23][24][25][26] , there are evidence of incoherency between the classical theory of photoconductance and experiments, such as the unexpected decrease of gain with light intensity 17,21,24 or its dependency with temperature [START_REF] Nie | [END_REF] . This paper proposes, based on Drift Diffusion theory, an original explanation of the origin of gain in vertical photoconductor. In particular, it demonstrates that light absorption can trigger charges injection by the contacts (a mechanism that is named in this paper "gain by injection enhancement") and discuss the conditions required to obtain this effect. It also explains why this phenomenon is much more pronounced at low light intensity, making the responsivity dependent of the light flux. In addition, the role of trapping as a booster of gain is analyzed. This theory is in very good agreement with new experiments performed on organic photoconductors. Moreover, the impact of gain on overall device performances (detectivity and cut-off frequency) is also discussed. Finally, several original solutions to enhance further the detectivity of vertical photoconductor are investigated. This work considers exclusively the case of organic photodetector with metallic contacts, even if we believe that it can be easily extended to other kind of materials or contacts, provided that specific effects, such as ionic transport in the case of perovskite material, or tunneling in the case of insulating contacts, are included. This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

2) Investigation of the main features of gain in organic photodetectors by numerical simulations

In this first section, the main features of gain in both organic photoconductors and photodiodes are investigated by Drift-Diffusion simulations. Details about devices, materials parameters and models used in the simulation are presented first, before discussion. The device structures of the organic photodetectors considered in this work are depicted in Figure 1. In all devices, the photoactive layer is a poly(2,7-carbazole-alt-4,7-dithienyl-2,1,3-benzothiadiazole) (PCDTBT) and phenyl-C61-butyric acid methyl ester (PC60BM) PCDTBT:PC 60 BM bulk heterojunction, materials known to lead to overall good device performances and well documented in the literature 8 . The blend was treated as a homogeneous semiconductor, following the approach proposed by Koster et al. [START_REF] Koster | [END_REF] . Three typical template device architectures were considered. The first one is a photoconductor with symmetrical hole injecting contacts with electrode work function at c.a. -5.4 eV, typically the poly(3,4-ethylenedioxythiophene) polystyrene sulfonate (PEDOT:PSS). This architecture leads to a "hole-only" or "p-only" device. For the second architecture, mid gap symmetrical contacts with electrode work function at c.a. -4.4 eV (typically Aluminum or ITO with PEIE (Polyethylenimine ethoxylated 32 )) have been used. The third architecture is a photodiode with asymmetric contacts. It combines each of the previously mentioned contacts for comparison, i.e. one hole and one electron injecting contact (even though a better electron injecting contact would have been used in typical photodiode, however without affecting the generality of the results). In analogy with doped crystalline photoconductors, the "p-only" device is called extrinsic photoconductor, while the device with mid gap contact is called intrinsic photoconductor (even if in our case, no doping was used).

The transport of injected and photo-generated carriers within the device was modelled through the Drift-Diffusion approximation, which has been shown by several authors 9,[START_REF] Koster | [END_REF][33][34][35] to model efficiently either injection, transport, generation or trapping of carriers both in the dark and illuminated regimes. The numerical simulations were performed using the Fluxim software. [START_REF] Fluxim | SETFOS: Semiconducting emissive thin film optics simulator software[END_REF] It allows the introduction of experimentally determined parameters e.g. electrode work function, charge carrier mobility, etc. This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

Importantly, the classical boundary conditions introduced by Schottky [START_REF] Sze | Physics of Semiconductor Devices[END_REF] were used, which consists in keeping the electron and hole concentrations at electrode contacts constant under field and lighting and equal to their equilibrium values. Previous results have shown that these boundary conditions were in good agreement with experiments, accurately predicting the impact of contact work function modification and temperature [START_REF] Kielar | [END_REF] . Moreover, the HOMO and LUMO effective density of states N 0 was selected according to a common value found in literature 38 to 10 26 m -3 . The blend material parameters and thickness (300 nm) are kept equal in the three devices and are reported in Table 1.

As in PCDTBT:PC 60 BM bulk heterojunction, the mobility of electrons and holes have been found of the same order of magnitude 39 and in order to simplify the discussion on the origin of gain, both electron and hole mobilities were also taken constant and equal to 5x10 -3 cm 2 /V.s. Of course, it is also possible to consider different mobility models for electrons and for holes and their temperature and electric field dependency. It will be done later on in this paper.

As it has been suggested in the literature that the gain was impacted by electron trap levels lower than the LUMO level (shallow trap), simulations were performed with and without shallow traps. In many organic photodiodes (including PCDTBT:PC 60 BM blends), these traps are typically observed in oxygen contaminated bulk heterojunction 37 40,41 . In absence of light, these traps are empty. In presence of light, photo generated electron may remain trapped within the traps, creating negative charges within the active layer of the device. The C n , C p trap parameters are not accurately known, and a large set of value can be found in the literature [START_REF] Kielar | [END_REF]42 . We decided arbitrary to use the value of 10 -13 cm 3 /s. However, note that the recombination lifetime ((C n +C p ) N T ) -1 considered in this work for trap concentration in the range of 10 17 -10 19 cm -3 are in good agreement with experimental measurements (50 µs -0.5µs) 42,43 .

All device considered are vertical devices, where the light is entering through one semitransparent contact (and not laterally in conventional crystalline photoconductors). To model the light propagation in the stack, the optical indexes of the PCDTBT:PC 60 BM were used. These values were extracted from previous reflectance and transmittance experiments, following the approach proposed in 44 . The monochromatic light source was set (430nm) with a constant power of 0.1 mW/ cm 2 . This wavelength was chosen to lead to a globally uniform generation in the organic blend of 300nm thickness in order to simplify the first step of the discussion. When comparing with experiments (section 4), a green LED was used.

Finally yet importantly, each photon absorbed in the blend is supposed to lead to the creation of a dissociated electron-hole pair, which, of course, is a very optimistic assumption, because of the difficulty of dissociating excitons in the real blend.

b) Results: Comparison between photodiodes and photoresistor

Main results are summarized in Table 2. As expected, the defect-free photodiode has the lower level of dark current, followed by the defect-free intrinsic photo-resistance and the photodiode with traps. Extrinsic photoresistor, with or without traps, have a larger dark current. Note that the value of the dark current are (as expected) largely dependent of the contact nature.

More surprising, only the extrinsic photo-resistance shows significant photo-multiplication (gain larger than unity). This gain is very important in presence of traps (∼1150), in agreement with previous experimental observation in the literature (see for instance the review 16 ). However, a gain is also present in absence of trap, even though at very modest level (∼2). Clearly, this low value of gain has certainly low interest for application.

However, the observation of gain in absence of traps is interesting from a fundamental point of view, as it proves that photomultiplication is not the direct consequence of trapping. In addition, it cannot be explained with the classical theory of gain [START_REF] Sze | Physics of Semiconductor Devices[END_REF] in conventional planar photo-resistance. Indeed, according to this model, the gain is equal to:

𝐺𝑎𝑖𝑛 𝑐𝑙𝑎𝑠𝑠 = 𝜏 𝜏 𝑐𝑜𝑙𝑙 (1)
Where  is the charge carrier lifetime due to recombination, and 𝜏 𝑐𝑜𝑙𝑙 the collection or transit time, given by:

𝜏 𝑐𝑜𝑙𝑙 = ( µ 𝑛 𝑉 𝑎 𝐿 2 + µ 𝑝 𝑉 𝑎 𝐿 2 ) -1 (2)
Where µ n (resp. µ p ) is the electron (resp. hole) mobility, L is the active layer thickness, V a the applied voltage. In absence of recombination ( tends to the infinite), the gain Gain class should diverge, suggesting that this theory does not apply in absence of recombination.

In conclusion, simulations confirm that photo-multiplication may occur extrinsic photo-resistance. This mechanism is particularly strong in presence of electrons traps (for hole only photo-resistance), but may exist also in principle in absence of any traps. This result suggests that traps only amplifies photo-multiplications, but are not necessarily required for it to occur. Table 2: Summary of the results of the simulations of the first section, for simplified PCDTBT:PC 60 BM , 300nm thick, photodiode, intrinsic and extrinsic photo-resistance considered at -2V with a constant illumination of 430nm set at 0.1mW/cm 2 , traps located at 0.7eV below the LUMO level with a density of 10 19 cm -3 , and a constant mobility of 0.005cm²/V.s (equal for electrons and holes). The complete simulation parameters are summarized in Table 1. c) Results: Impact of light irradiance on the gain of extrinsic photo-resistance.

Notation

In order to get more insight on the photo-multiplication effect, the impact of the irradiance level on the value of gain was investigated for the three device architectures mentioned above (Figure 2). Devices with no traps and no photo-multiplication (gain closed to unity), i.e. the trap free photodiodes and intrinsic photo-resistances, does not show, as expected, any irradiance dependency. In presence of traps, the responsivity of these devices is reduced, and further degraded when irradiance increases, until reaching a constant level. This effect is expected, as it is consistent with the progressive process of trap filling until trap saturation.

Interestingly, for device featuring photo-multiplication, the gain is constant at low level of irradiance and then decreases whenever traps are present or not (even if this droop is much less pronounced for the photo-resistance without traps). This trend is intriguing, as the gain decreases despite the fact that traps are more and more filled by photogeneration. It has already been observed in several previous experimental papers 17,21,23,24 .

In any case, a high value of responsivity, and its decrease with light irradiance can be seen as an indication of the presence of gain, in agreement with the theory discussed in this paper, allowing to discriminate from other potential source of gain. In particular, the classical theory of gain in conventional planar photo-resistance does not predict any similar irradiance level dependency.

Even if it is also based on a drift diffusion approach, we believe that the classical theory of gain in planar photo-resistance does not apply in vertical devices. In fact, this analytical theory explicitly assumes that the conductance in presence of light remains the same for all applied voltage, which may be true when the light does not enter through the electrodes, but by a lateral side. In our case however, simulations are not based on any a-priori approximation and leads, as discussed, to different results. The simulation was performed for a monochromatic illumination set at 430nm. The illumination intensity power varies from 10 -4 W/m 2 to 10 4 W/m 2 . The other simulation parameters are kept identical to those used in the previous section and are reported in Table 1.

The gain is also voltage and wavelength dependent, as discussed in the supplementary materials. This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

3) Theory of gain by injection enhancement: model and discussion

As seen previously, steady state drift diffusion simulations confirm the existence of photomultiplication, however it is not obvious just from simulation to get detailed insight on the origin of the phenomenon. Indeed, as drift and diffusion equations compute the total electron and hole concentration at each point, it is not possible to discriminate from the direct result of calculations carriers that have been directly photo-generated by photon absorption from other carriers induced by photo-multiplication.

For this reason, a model is proposed in this section, allowing to investigate the origin of photomultiplication. Precious indications can be obtained from this model, as seen later. According to the calculations detailed in appendix A, it is possible to get an expression of the gain (noted "Gain" to avoid confusion with the generation rate G) in the framework of the drift diffusion model, without any approximations. This expression is:

𝐺𝑎𝑖𝑛 = | < 𝜙 𝐺 (𝑥) -𝜙 𝑅 (𝑥) > 𝑖 -< 𝜙 𝑅𝑑 (𝑥) > 𝑑 𝜙 𝐺𝑡 + -𝛽 𝑛 𝜙 ̃𝑛𝑑 + 𝛽 𝑝 𝜙 ̃𝑝𝑑 𝜙 𝐺𝑡 | (3)
where 𝜙 ̃𝑛𝑑 (resp. 𝜙 ̃𝑝𝑑 ) is the electron (resp. hole) carrier flux in dark condition, 𝜙 𝐺𝑡 is the photon flux (m -2 s -1 ) absorbed in the active layer (thickness L), and

𝜙 𝐺 (𝑥) = ∫ 𝐺(𝑥 ′ )𝑑𝑥′ 𝑥 0 (4)
𝜙 𝑅 (𝑥) = ∫ 𝑅(𝑥 ′ )𝑑𝑥′ 𝑥 0 (5)
𝜙 𝐺𝑡 = 𝜙 𝐺 (𝐿) = ∫ 𝐺(𝑥 ′ )𝑑𝑥′ 𝐿 0 (6)
G (resp. R) is the generation (resp. recombination) rate under illumination. See the appendix for the definition of other parameters. This formula allows to recalculate the gain using the electrostatical potential V and V d extracted from simulation in very good agreement with direct simulation results (see Table 3), and more importantly to identify its different contributions, as discussed in the following. However, let us stress that equation (3) cannot allow to calculate the gain alone, as the potential profile V under illumination and V d under dark are not known, and have to be extracted from simulations. The aim of this model is not to calculate the gain, but rather to help understanding simulation results, as detailed in the following. A complete analytical model of gain would require formula for V(x) and V d (x), but this is beyond the scope of this paper. This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. 3) and compared with data directly extracted from simulation) for a PCDTBT:PC 60 BM photodiode, photo-resistance intrinsic and photo-resistance extrinsic at -2V. The other simulation parameters are kept identical to those used in the previous section reported in Table 2.

PLEASE CITE THIS ARTICLE AS

Let us consider the first term, , defined as:

 = |< 𝜙 𝐺 (𝑥) -𝜙 𝑅 (𝑥) > 𝑖 -< 𝜙 𝑅𝑑 (𝑥) > 𝑑 | 𝜙 𝐺 (𝐿) (9)
This  term would correspond to the internal quantum efficiency in absence of photo-multiplication. Indeed, let us prove that  cannot be higher than 1. As recombination tends to decrease , we have:

 < |< 𝜙 𝐺 (𝑥) > 𝑖 | 𝜙 𝐺 (𝐿) < max(𝐴, 𝐵) (10)
Moreover, both A and B are lower than 1, as 𝜙 𝐺 (𝑥) < 𝜙 𝐺 (𝐿). Being lower than unity, this term cannot induce gain.

The second term, under some circumstances, may be responsible for photo-multiplication, as explained in the following. First of all, let us stress that this term is proportional to the dark current, as 𝜙 ̃𝑛𝑑 (resp. 𝜙 ̃𝑝𝑑 ) is the flux of electrons under dark condition (resp. the flux of holes). First interesting conclusion, the photo-multiplication process only occurs in device featuring significant dark current, which explains why photodiodes or intrinsic photo-resistance have not shown any gain in previous simulations (see Table 3). Moreover, being proportional to the dark electron or holes flux, this term represents clearly an additional source of current provided by contact injection. As discussed in the following, this condition is necessary, but not sufficient. This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

As extrinsic photo-resistances are either P or N only devices (hole or electron are majority carrier), there is no need to consider both terms simultaneously. Let us consider the case of P only device. In this case, the gain equation simplifies:

𝐺𝑎𝑖𝑛 = |𝛼 + 𝛽 𝑝 𝜙 ̃𝑝𝑑 𝜙 𝐺𝑡 | (13)
The  p coefficient is given by:

𝛽 𝑝 = ∫ exp ( 𝑒𝑉 𝑑 (𝑥) 𝑘𝑇 ) 𝑑𝑥 𝐿 0 ∫ exp ( 𝑒𝑉(𝑥) 𝑘𝑇 ) 𝑑𝑥 𝐿 0 -1 = ∫ exp ( 𝑒𝑉 𝑑 (𝑥) 𝑘𝑇 ) 𝑑𝑥 𝐿 0 ∫ exp ( 𝑒𝑉 𝑑 (𝑥) 𝑘𝑇 ) exp ( 𝑒Δ𝑉(𝑥) 𝑘𝑇 ) 𝑑𝑥 𝐿 0 -1 (14)
where

∆V(x) = V(x) -V d (x)
is the difference of electrostatic potential in presence of light (V) or in absence of light (V d ). A simple mathematical analysis of equation ( 14) reveals that 𝛽 𝑝 is necessarily higher than -1, and can take large positive value when ∆V << -kT. This later case corresponds to device with large gain. However, this condition is not easy to achieve, as discussed in the following. ∆

V ≠ 0 means that light carrier generation should induce a significant modification of the charge within the device. This condition requires particular conditions to be satisfied. Indeed, it is generally assume on the contrary that there is no significant change in the total charge [START_REF] Bube | Photoconductivity of Solids[END_REF] after illumination ("quasi neutrality" approximation), as both negative charges (electrons) and positive charges (holes) are generated simultaneously. The Photo-Dember effect [START_REF] Sze | Physics of Semiconductor Devices[END_REF] is an example of known situation where ∆V ≠ 0 on the contrary, due to the difference of mobility of photogenerated electrons and holes. In our case, if the difference of mobility may also indeed induce gain, the most efficient mechanism remains trapping. Indeed, traps can capture a significant amount of one type of carrier, inducing a charge difference between dark and light conditions, explaining why traps are so efficient to create gain.

A closer analysis reveals that additional requirements are needed for the charge modification induced by photogenerated carriers to create significant gain. Indeed, first ∆V need to be negative, then ∆V(x) need to take value at a suitable position (depending on V d (x)), in order to contribute significantly to the integral of Eq. ( 14). These unexpected requirements can be understood as follow: as the photomultiplicated carriers came from the contacts, the additional electrical field induced by light absorption need to reinforce the initial dark electrical field (and not decrease it), in order to trigger injection (at least in one region of the active layer). This situation is illustrated in Figure 3. In this figure, the excess of dark electric field E d (x, -2V) -E d (x, 0V) is plotted together with the additional field induced by light E(x, -2V) -E d (x, -2V). The excess of dark field is responsible for hole injection, causing dark current 𝑒𝜙p 𝑑 . It is largely negative, meaning that holes are moving from right (x=L) to left (x=0). This device features a large level of gain because light induces (by electron trapping) a large spike of negative field close to the right electrode (x=L), reinforcing injection. In consequence, the spatial repartition of the generation rate (which has been considered relatively uniform so far) is expected to have an impact on the gain. An example is shown in the supplementary material.

Trap parameters used in Figure 3 E 4.

All these explanations conclude that the photo-multiplication process is the consequence of an additional injection of majority carriers induced by electrostatic. We call it "gain by injection enhancement".

The dependency of gain with light flux (reported in Figure 2 for instance) can be understood as follow. First of all, let us note that the two terms in equation ( 13) depend on 𝜙 𝐺𝑡 , the total number of generated electron-hole pair: the ratio 𝜙 ̃𝑝𝑑 /𝜙 𝐺𝑡 and 𝛽 𝑝 . If 𝛽 𝑝 linearly increase with 𝜙 𝐺𝑡 , the gain is constant and if 𝛽 𝑝 saturates the gain decrease as 1/𝜙 𝐺𝑡 . Keeping this in mind, the trend of gain versus irradiance can be explained as follow when trapping is involved. At low irradiance flux, traps are mostly empty and keep on charging when irradiance increases. This results on a linear increase of 𝛽 𝑝 with 𝜙 𝐺𝑡 , leading to a constant gain. At high irradiance however, the trap occupancy saturates, resulting in a constant 𝛽 𝑝 and thus a decreasing gain.

4) Experimental validation

In order to validate our theory by comparison with experiments, PCDTBT:PC 60 BM photoresistor have been processed, using a process very similar to the one presented in reference 8 . For this study, PEDOT:PSS have been used for the top (~1.5 µm) and bottom electrodes (~50 nm) to have symmetric architecture. The active layer thickness was measured approximately equal to 300 nm. For a comparison between our model and experimental data, the following field dependent mobility model has been used (assuming same mobility for electrons and holes in agreement with the literature 39 ).

µ(𝐸) = µ 𝑝 exp (𝐵√𝐸) (15) Interestingly, all parameters needed for reproducing experiments were already determined in a previous work [START_REF] Kielar | [END_REF] , by comparison between simulations and experiments on photodiodes (not photoresistor). No further adjustment was needed to reproduce photoresistor results, except a minor adjustment of the value of the work functions of PEDOT electrodes (around 5.36 eV). The bottom PEDOT:PSS electrode is spin-coated whereas the top PEDOT:PSS electrode is screen-printed. Moreover, the deep traps present in the simulation, and responsible for gain in this case, were natively This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

present in the active layer, and already identified in reference [START_REF] Kielar | [END_REF] , in order to reproduce correctly the level on illumination current on virgin and aged devices. The same trap parameters were used for both photodiode and photoresistor simulations. In experiments, the light source was a green LED (536 nm, close to the maximum sensitivity peak of PCDTBT:PC 60 BM photodiode 8 ). The gain in these devices was found equal to ∼5 at -2V @ 24.671 W/m² (by simulations, see Figure 4 (b)). 5.

Simulations parameters

The agreement between experiments and simulation is satisfactory, in I-V experiments (Figure 4) as well as in responsivity versus irradiance (Figure 5), which confirms the validity on the simulation and theory discussed in the paper. In particular, the decrease of responsivity with irradiance agrees with the theoretical results presented in Figure 2. This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

5) Consequences on device performances

In the previous section, the physical mechanism of gain in photoresistor by injection enhancement has been investigated. In this section, the impact of this gain on device performance is discussed and compared with other devices.

a) Detectivity

One of the main figure of merit of photodetectors is detectivity 3 D*. A device featuring a high detectivity has typically a good responsivity and a low noise level, allowing detecting low irradiance signal. The high gain at low irradiance of photoresistor enhanced by photomultiplication is clearly an advantage for the detectivity. However, it comes with a price: as previously discussed, a high gain implies a high level of dark current, leading to a high level of dark current noise. This fact rises the following question: is there a trade-off between gain and dark current noise, allowing photoresistor to be competitive with photodiodes?

To answer this question, the detectivity of several devices has been plotted versus trap concentrations on Figure 6. Three devices with no photo-multiplication (two photodiodes and an intrinsic photoresistor) were compared with two extrinsic photoresistors presenting photo-multiplication due to trapping. Devices are essentially similar to those reported in Table 1. All devices are composed of a blend of PCDTBT:PC 60 BM, and the traps used in the simulations were the same (volume trap uniformly distributed, electron acceptor traps, mid gap, with Cn = Cp = 10 -13 (cm 3 /s)).

The two photodiodes considered have the same electrodes: a hole injecting contact -5.1 eV and an electron injecting contact -4 eV. However, the first one (named PD) is assumed ideal. It means that the dark current due to contact injection is extremely low. Recombination does not increase it sufficiently, and for this reason, the detectivity remains photon shot noise limited, independent of the trap concentration, extremely high and unrealistic (10 15 Jones). The best reported organic photodiodes, as in ref 8 for instance have a more reasonable detectivity in the range 3.10 13 Jones, due to an enhanced dark current induced by internal shunt. For this reason, a shunt has been included in our simulation to reproduce state of the art photodiodes (PD(Rp) symbols). In this case, the detectivity decreases with the concentration of recombination traps, essentially because the responsivity is degraded by trapping (shunt is assumed independent of trapping).

The intrinsic photoresistor has mid gap electrodes, leading to an extremely low dark current in absence of traps. In presence of traps however, the dark current increases and the responsivity decreases, leading, as expected, to a significant degradation of the detectivity.

Finally, two extrinsic P-only photoresistors were also considered (work function -5.1 eV and -5.3 eV).

Being single carrier devices, there is no significant evolution of the dark current by recombination (no electrons in dark). However, the increase of trap concentration induces a significant increase of the gain, explaining why the detectivity improves with the trap concentration in both case. The slightly better detectivity of the 5.1 eV photoresistor is explained by its better value of  h . This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. The main conclusion of these simulations is that the detectivity of photoresistors can reach good value (> 10 12 Jones) despite their high level of dark current. However, these theoretical performances remains lower that the best value reported in the case of state of the art visible organic photodiodes (3.10 13 Jones). In consequence, photoresistor enhanced by gain may be an interesting alternative in specific technology where it is difficult to achieve good photodiodes, either because contacts electrodes are poor, or because it is difficult to reduce traps concentration (due to impurities and disorder for instance). Also, photoresistances might be more resilient to contamination than photodiodes, as contamination might increase the traps concentration responsible for gain. However, this point has not been investigated in this work and would require further attention.

PLEASE CITE THIS ARTICLE AS

b) Frequency operation

Another figure of merit of photodetector is its bandwidth. In order to investigate the frequency operation of photo-resistance featuring photo-multiplication due to traps, simulations were performed, fixing the applied voltage, and modulating the frequency dependency of the illumination. The device considered is the same as in the experimental section (parameters listed in Table 5), including field dependent mobility.

In Figure 7, we can observe two frequency cut-offs. At about 1 Hz, the gain rapidly drops off. It corresponds to the slow response of deep traps. Between 1 Hz and 1 MHz, the device has a low gain value, which approximately correspond to its gain in absence of traps. Above 1 MHz, the device has no response. It corresponds to its time of flight cut-off frequency (mobility is enhanced in this case by the large applied voltage and intrinsic RC time is quite short). As expected, the traps concentration affects the value of the steady state gain, but not its frequency operation.

Other simulations indicate that the gain frequency cut off depend on the trap level position (shallow traps are faster) but in any case, the first cut off frequency never exceed hundreds of Hz, in agreement with experimental results reviewed in details in 3 .

Our simulations thus confirm that the gain by injection enhancement due to traps is a rather slow process, penalizing potential application for this type of devices. This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

Figure 7 : Gain dependency versus frequency (applied voltage -2V) in a P only extrinsic photoresistor with different traps concentration (Fig. 7a) and trap position (Fig. 7b). Parameters are identical to those reported in Table 2.

6) Conclusions

The mechanism of Trap Assisted Photo-Multiplication in vertical organic photoresistor, observed in several previous experimental works, was investigated on a theoretical point of view in this paper, providing several interesting conclusions.

First of all, drift diffusion simulations have confirmed that a gain may occur, essentially in P only (or N only) photoresistor, but not in classical photodiodes or intrinsic photoresistor. Our model has allowed to analyze in details the origin of this gain. It turns out that photo-generated carriers, if able to induce a significant modification of electric field, may trigger additional carrier injection from the contacts. This mechanism has been called "gain by injection enhancement". It is responsible for a significant dependency of the photoresistor responsivity with irradiance (where "good" photodiodes have typically a constant responsivity), this effect being more sensitive at low irradiance. If carrier injection is limited, as in photodiodes or in intrinsic photoresistor, the gain remains negligible, but can exceed unity in the case of extrinsic photoresistor. Traps, fixing photogenerated carriers in the active layer, (electrons for P only devices, or hole for N only devices) may be an excellent way to induce such additional electric field, even if this mechanism may appear also in absence of traps in principle (however at moderate level). If the traps considered in this work were uniformly distributed in the active layer, it is clear however that further optimization of this mechanism is possible by playing with the localization and the nature of traps.

All these conclusions were confirmed by comparison with experimental P only PCDTBT:PC60BM photoresistors. Interestingly, once material parameters (mobility, traps …) were extracted by fitting photodiodes experiments, we found that no further adjustment were required to reproduce satisfactorily photoresistor experiments, confirming the predictive nature of our simulations.

In term of device performances however, gain by injection enhancement has several limitations. First, when the gain is achieved by trapping, it is an extremely slow mechanism (depending on the trap energy, limited a hundred of Hz in the best case). Second, as gain takes interesting value only when dark currents are significant, it cannot achieve very high value of specific detectivity (it cannot reach the photon shot noise limit for instance). In the template material considered in our work, (PCDTBT:PC60BM blends) photoresistor simulations indicate that it was not possible to reach the value of detectivity already obtained in the best experimental results reported using photodiodes. However, it may remain attractive for other technology and material, where the nature of material and contacts does not allow to scale down the dark current of photodiodes. This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

Note that  G (L) is the total number of electrons holes pairs generated, thus the number of photons absorbed.

In dark condition, we have :

𝜙 𝑛𝑑 (𝑥) = 𝜙 ̃𝑛𝑑 -𝜙 𝑅𝑑 (𝑥) (

  a) Description of the template devices and model

Figure 1 :

 1 Figure 1: a) Schematic view of the organic photo-resistance considered in this section. b) Energy diagram of the PCDTBT:PC 60 BM photo-resistance at equilibrium and zero bias

  PD No trap PR(ext) With trap PR(ext) No trap PR(int) With trap PR(int) No trap

Figure 2 :

 2 Figure2: Gain versus the irradiance for a PCDTBT:PC 60 BM different structure at -2V. The simulation was performed for a monochromatic illumination set at 430nm. The illumination intensity power varies from 10 -4 W/m 2 to 10 4 W/m 2 . The other simulation parameters are kept identical to those used in the previous section and are reported in Table1.

  𝐺 (𝑥)/𝜙 𝐺 (𝐿)] exp (

Figure 3 :

 3 Figure 3: Simulated electric field difference between illumination and dark condition (E(x) -E d (x)) (blue circle). The black triangle curve corresponds to the electric field variation between -2V and 0V in dark condition only. The simulation parameters are reported in Table4.

Figure 4 :

 4 Figure 4: Illuminated and dark I-V curves, dots are experimental data, solid lines are simulations. (a) Comparison between simulation and experiments using the parameter of Table 5. (b) Corresponding calculated gain versus the voltage.

Figure 5 :

 5 Figure 5: Responsivity versus irradiance for two different applied voltages (a) V a = -2V and (b) V a = -0.5V. Dots are experimental data, solid lines are simulations. The simulations were performed with the parameters listed in Table5.

7 Figure 6 :

 76 Figure6: Detectivity versus trap density for a PCDTBT:PC 60 BM different devices at 0.5V. "PD" refers to an ideal photodiode, "PD(Rp)" refers to a photodiode where the dark current is limited by a shunt resistance (see the text for details), "PR(ext)" refers to extrinsic photoresistor with gain due to traps, "PR(int)" refers to a instrinsic photoresistor (work function at 4.4 eV) with no gain.

<< 1 0 10 10 10 10 10 7 10

 11010101010710 A13) This equation can be compacted introducing the following notation: 𝜙 ̃𝑛 -𝜙 ̃𝑛𝑑 = 𝛽 𝑛 𝜙 ̃𝑛𝑑 -< 𝜙 𝐺 (𝑥) -𝜙 𝑅 (𝑥) > 𝑛𝑖 -< 𝜙 𝑅𝑑 (𝑥) > 𝑛𝑑 𝜙 𝐺 (𝑥) -𝜙 𝑅 (𝑥) > 𝑛𝑖 = ∫ [𝜙 𝐺 (𝑥) -𝜙 𝑅 (𝑥)] exp (procedure for holes conservation equation, we found : 𝜙 ̃𝑝 -𝜙 ̃𝑝𝑑 = 𝛽 𝑝 𝜙 ̃𝑝𝑑 -< 𝜙 𝐺 (𝑥) -𝜙 𝑅 (𝑥) > 𝑝𝑖 -< 𝜙 𝑅𝑑 (𝑥) > 𝑝𝑑 𝜙 𝐺 (𝑥) -𝜙 𝑅 (𝑥) > 𝑝𝑖 = ∫ [𝜙 𝐺 (𝑥) -𝜙 𝑅 (𝑥)] exp ( 𝑒𝑉This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1Electric field variation under illumination (kV/cm) Electric Field Variation under illumination (kV/cm)This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has bePLEASE CITE THIS ARTICLE AS DOI:This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has bePLEASE CITE THIS ARTICLE AS DOI:This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has bePLEASE CITE THIS ARTICLE AS DOI:This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has bePLEASE CITE THIS ARTICLE AS DOI:This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has bePLEASE CITE THIS ARTICLE AS DOI:This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has bePLEASE CITE THIS ARTICLE AS DOI:This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has be PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

Table 1 :

 1 This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. Summary of the simulation parameters used in this section.

	Global parameters
	HOMO	3.7eV
	LUMO	5.5eV
	N 0	10 26 cm -3
	Mobility (µ n =µ p )	5x10 -3 cm 2 /V.s
	Work function	Between 4.4eV and 5.36eV
	Wavelength	430nm (corresponding to an absorbance of 0.842)
	Intensity power	0.1mW/cm²
	Length/thickness	300nm
	Trap parameters (electron acceptor)
	E T	0.7 below the LUMO level
	N T (N type)	1x10 19 cm -3
	Capture rate: Cp, Cn	1x10 -13 cm 3 /s

PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

Dark current (mA/cm²) Electron current (mA/cm²) Hole current (mA/cm²) S (A/W) Gain

	Photodiode	PD	10 -11	0	10 -11	0.287	0.982
	Photodiode with trap	PD(T)	10 -10	0	10 -10	0.117	0.4
	intrinsic photo-resistance no trap	PR(int)	9.2x10 -6	9.2x10 -6	1.8x10 -12	0.282	0.966
	intrinsic photo-resistance with trap	PR(int)T	10 -10	0	10 -10	0.06	0.2
	extrinsic photo-resistance no trap	PR(ext)	318	0	214	0.6	1.93
	extrinsic photo-resistance with trap	PR(ext)T	318	0	318	336	1150

DOI: 10.1063/1.5126338 (mA/cm²)

	Photodiode	0.98 -4x10 -6	0	3x10 -5	-10 -11	0 10 -11	0.98	0.98
	Photodiode with trap	0.41	-0.451	0	0.8	-10 10	0 5x10 -9	0.4	0.41
	photo-resistance intrinsic (No trap)	0.97	-10 -5	9x10 -6	6x10 -5	10 -12	3x10 -4 6x10 -11	0.966	0.97
	photo-resistance extrinsic (No trap)	0.93	4x10 -7	0	2x10 -5	214	0 5x10 4	1.93	1.93
	photo-resistance intrinsic (Trap)	0.20	5x10 -5	0	-1	5x10 -9	0 5x10 -9	0.2	0.20
	photo-resistance extrinsic (Trap)	0.92	10 -3	0	0.105	318	0 10 4	1150	1156

Table 3 :

 3 Summary of the result of the simulation and coefficient β e , β h and α (The value and gain was recalculated with formula (

Table 4 :

 4 This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset. Summary of the trap parameters used in Figure3.

	T	0.7 below the LUMO level
		(electron acceptor)
	Capture rate: Cp, Cn	1x10 -13 cm 3 /s
	N T	110 19 (cm -3 )
	β h	0.105
	β e	1.3x10 -3
	Gain	1150

PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

Table 5 :

 5 Parameters used in simulations to reproduce photoresistor experiments

	Work function	5.36 eV	5.36 eV
	µ 0	µ n =µ p =810 -3 cm²/V.s
	B	0.6510 -2 (cm/V) 0.5
		536 nm
	Irradiance	24.671 (W/m²)
	Active layer thickness	300 nm
	Active area	2.510 -6 (m 2 )
	Trap parameters
	E T (eV)	0.7eV below the LUMO level
	N T	310 17 (cm -3 )
	C p /C

n 10 -13 (cm 3 /s)

This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has be PLEASE CITE THIS ARTICLE AS DOI: 10.1063/1.5126338

7) Supplementary material

More information about the impact of applied voltage and wavelength is given in the supplementary material.

8) Acknowledgment

This work has been funded by the French National Research Agency (ANR) through the TAPIR project N° ANR-15-CE24-0024-01.

Appendix A : Derivation of the formula of gain by injection enhancement

According to the drift diffusion model, the electron and hole carrier flux (m -2 s -1 ) at each point x are given by:

where µ n and µ p are the electrons and holes mobility, k the Boltzmann constant, T the temperature, n and p the carrier concentration (electrons and holes respectively), e the elementary electron charge and V the potential.

The conservation equation of electrons and holes leads to :

where G(x) and R(x) are respectively the electron-hole carrier generation and recombination rates. These equations ensure the conservation of the total current density j illum (Am -2 ), as :

The conservation equation of electron can be written [START_REF] Sze | Physics of Semiconductor Devices[END_REF] in the following form by multiplication of each side par exp(-eV/kT):

By integrating this equation through the entire structure and using the classical Schottky boundary conditions (the same boundary condition implemented in the simulation), the right-hand side simplifies and the conservation equation becomes:

where L is the active layer thicknesses, n L = n(x=L) (resp. n 0 = n(x=0)) is the electron concentration at the left (resp. right) contact, and V L -V 0 = V a , the applied voltage.

Therefore, thanks to the boundary conditions, this integral has the same value in dark and illumination conditions. Consequently,

where V d (x) is the potential profile in dark condition, φ nd (x) the flux of electrons in dark.

According Eq.(A3), φ n (x) can be replaced by: According Eq. (A4), the constant current density under illumination is given by : 𝑗 𝑖𝑙𝑙𝑢𝑚 = -𝑒𝜙 𝑛 (𝑥) + 𝑒𝜙 𝑝 (𝑥) = -𝑒𝜙 𝑛 (0) + 𝑒𝜙 𝑝 (0) = -𝑒𝜙 ̃𝑛 + 𝑒𝜙 ̃𝑝 (A22) Similarly, in dark condition: 𝑗 𝑑𝑎𝑟𝑘 = -𝑒𝜙 𝑛𝑑 (𝑥) + 𝑒𝜙 𝑝𝑑 (𝑥) = -𝑒𝜙 𝑛𝑑 (0) + 𝑒𝜙 𝑝𝑑 (0) = -𝑒𝜙 ̃𝑛𝑑 + 𝑒𝜙 ̃𝑝𝑑 (A23) Thus, the photocurrent 𝑗 𝑖𝑙𝑙𝑢𝑚 -𝑗 𝑑𝑎𝑟𝑘 is given by :

The previous equation can compacted as :

Finally, the Gain is defined as the number of photo electrons collected, divided by the number of photon absorbed:

Using Eq. (A25), we obtain the following equation for the Gain :

This is the author's peer reviewed, accepted manuscript. However, the online version of record will be different from this version once it has been copyedited and typeset.

PLEASE CITE THIS ARTICLE AS