Here K is the constant of integration. Solving the quadratic equation for ˙ p , we obtain two first order ODEs, The rate of change of the internal energy, ˙ E , is the sum of two parts, namely, a work flux ˙ W and a heat flux ˙ Q , Note that the particle exchange contributes to the heat flux [last term in Eq. (6b)]. When the energy level is below the chemical potential, the direction of heat flow is opposite to the direction of tunnelling. The net total work and net total heat during the pro- cess of duration τ are obtained as functionals of the oc- cupation probability, The upper sign ( − ) should be used for upward processes in which the quantum level is raised and the lower sign (+) for downward processes. It is worth mentioning a useful symmetry between electrons and holes. We are using the state of an electron, ϵ ( t ) and p ( t ), to describe the state of system. Instead, we can also use the state of holes, − ϵ ( t ) and 1 − p ( t ). If p ( t ) is a solution for an upward process, then 1 p ( t ) is a solution for a downward process with − ϵ ( t ). Hence we do not need to calculate the downward process separately, as it follows from this symmetry. For cyclic processes, we have E (0) = E ( τ ) and hence W + Q = 0. For mathematical simplicity, we evaluate power using net heat instead of net work: Before turning to the solution of the differential equa- tion (12), we examine the physical meaning of the con- stant K . Eliminating ˙ p in Eq. (11) by using the master equation (3), the resulting quadratic equation for p ( t ) leads to the relation Our goal is to maximize the power output and to eval- uate the corresponding efficiency. Power is a complicated functional of the time-dependent protocols in stage I and III, and an exact analytical analysis looks difficult at first sight. The optimization can however be done in two steps. First, we fix parameter values, τ$_{c}$ , τ$_{h}$ , p$_{0}$ and p$_{1}$ and maximize the power with respect to the functional space of ϵ ( t ). Since the total operation time is fixed, we just need to maximize the heat. Next, we further maxi- mize the power with respect to the remaining degrees of freedom τ$_{c}$ , τ$_{h}$ , p$_{0}$ and p$_{1}$ . The problem of maximizing heat or minimizing work for a single-level quantum dot moving between given initial and final energy states has already been analysed in detail in [14]. We reproduce the crucial steps of this analysis for self-consistency. This equation indicates that when K = 0, p ( t ) is the equi- librium distribution associated with the instantaneous value of the energy, implying that K = 0 corresponds to the quasi-static limit ( τ → ∞ ). As K increases, p ( t ) deviates from the equilibrium distribution. We conclude that K measures how far the state of the system deviates from the quasi static limit. We will use this insight below to obtain a perturbative solution for small dissipation by assuming that K is small. Next we proceed to solve Eq. (12). Separation of the variables p and t leads to the following explicit result for the upward processes: where To find the protocol that maximizes the heat, we do not search directly for the optimal schedule ϵ ( t ), but identify the optimal occupation probability p ( t ). This is done by expressing ϵ ( t ) in terms of p ( t ) and ˙ p ( t ), and rewriting the heat, Eq. (7b), as a functional of p ( t ) and ˙ p ( t ): Note that K depends solely on the operation time τ , the probabilities p$_{0}$ and p$_{1}$ , and the tunneling rate C but not on temperature. Unfortunately, the function (15) is quite complicated so we can not obtain an analytical expression The extremum is found via the standard Euler-Lagrange method, leading, after integration, to where β Q [ p ( · )] = ∫ τ 0 L ( p, ˙ p )d t, (9) L ≡ ln [ 1 p ( t ) + ˙ p ( t ) − 1 ] ˙ p ( t ) . (10) L − ˙ p ∂ L ∂ ˙ p = ˙ p 2 ( Cp + ˙ p )[ C (1 − p ) − ˙ p ] = K. (11) Cτ = F [ p ( τ ); K ] − F [ p (0); K ] . (16) F ( p ; K ) = − 1 2 ln p + 1 √$_{K}$ arctan [ 1 − 2 p √ K + 4 p (1 − p ) ] + 1 2 ln [ 2 p + K + √ K 2 + 4 Kp (1 − p ) 2(1 − p ) + K + √ K 2 + 4 Kp (1 − p ) ] . (15) Ct = F [ p ( t ); K ] − F [ p (0); K ] , (14) p ( t ) = 1 e βϵ ( t ) + 1 [ 1 + e βϵ ( t ) / 2 $^{√$_{K}$}$] . (13) P = − W τ = Q τ . (8) W [ p ( · )] = ∫ τ 0 ˙ ϵ ( t ) p ( t )d t (7a) III. OPTIMIZATION: GENERAL CASE ˙ p = K (1 − 2 p ) ∓ √ K 2 + 4 Kp (1 − p ) 2 C (1 + K ) . (12) ˙ W ≡ ˙ ϵp = ˙ εp − ˙ µp (6a) 3