Development of a High-Precision and Flexible Model for Accurate Simulation of Trolleybus Grids

The reduction of climate-changing emissions is vital, particularly in urban areas. To achieve this goal, the decarbonization of the transport sector is decisive. Trolleybus systems represent a solution towards sustainability as long as the integration of the electrical infrastructure with renewable sources as well as the introduction of full-electric fleets are actualized. To ensure a smart electrical transition the circuit modelling phase is crucial and must enable the most accurate simulations possible. The modularity of the catenary model and a versatile graphical user interface, which allows extensive topological change flexibility, distinguish the novel trolleybus network simulator presented in this study. The model completes the gaps left by earlier block-based simulation tools with its high precision and low processing effort. The suggested model’s enhancement in precision is evidenced by a graphical study of the voltage distribution examined in a Section of Bologna’s trolleybus grid. Finally, albeit non-iterative and thus computationally even cheaper, the proposed model is applicable for the simulation of more complex realistic scenarios involving the non-linearity of the equivalent power grid. The article will show satisfactory results on the convergence of the simulator.


I. INTRODUCTION
On the path towards smart city development, trolleybus systems are re-emerging as a springboard for carbon-neutral urban transport. In line with the global, green-oriented goals through 2050, public transport operators are being pushed to make existing infrastructures for electric traction sustainable and smart [1]. The main elements of a typical trolleybus electrical infrastructure are schematically depicted in Fig. 1. Trolleybus systems principally consist of DC traction power substations (TPSSs), e.g., uncontrolled rectifier substations (equipped with a twelve-pulse rectifier), and overhead con-The associate editor coordinating the review of this manuscript and approving it for publication was Arash Asrari . tact line (OCL) systems, or catenary. Differently from other electrified transport systems relying on dynamic conductive power transfer (e.g., tramways, a section of which is schematically illustrated in Fig. 1 for the sake of completeness, and metro systems), the power is carried by a two-wire OCL, so that the current return path is overhead. To guarantee the regularity of the traction system operation, OCLs are divided into feeding sections (FSs) that are separately powered in the case of fault [2]. Ordinary FSs are supplied bilaterally. In the event of bidirectional traffic, the two physically parallel OCLs are connected at different duly distributed points by equipotential bonding (EB), also known as voltage equalizers or current stabilizers, forming a double contact line (double-bifilar). The equipotentiality achieved by these electrical connections tends to limit the voltage drops along the catenary by reducing the electrical resistance, especially farther from the TPSSs.
Several actions can be identified to achieve greenhouse gas emission reduction targets. As proposed in the Trolley:2.0 project [3], the introduction of in-motion-charging (IMC) trolleybuses, for example, is a paramount step towards the independence from internal combustion engines. IMC trolleybuses are indeed vehicles powered by the OCL and equipped with battery storage systems to extend the journey in autonomous mode. The actions to mitigate climate-changing emissions can also be indirect, such as the introduction of higher-frequency cyclic bus schedules, which satisfies commuters' demand and hence encourages the use of public transport. However, an increase in service frequency as well as the operation of an IMC fleet, which requires more electric power than a conventional trolleybus fleet to feed the on-board batteries for the route extension, may excessively burden the traction network, leading to service disruptions in the worst-case scenario. Therefore, the whole infrastructure needs electrical support. Stationary energy storage systems (SESSs) may be considered, as shown in [4], along with renewable energy sources, to prevent undesired OCL voltage drops and to make the traction grid more eco-friendly.
Further clean air-oriented action concerns the deployment of electric vehicles (EVs). EVs are expected to be one of the leading players in the future green shift [5]. With a view to decreasing the charging times of EVs, the power of charging stations needs to be increased to 150 ÷ 350 kW. Currently, these powers can only be handled by DC fast charging stations (Level 3 chargers). Hence, catenary systems, being relatively powerful networks working in DC, may appear as a great opportunity to provide feeding points for interconnecting such state-of-the-art EV chargers [6], [7], [8], [9].
Monitoring the modernized trolleybus network through the technical innovations mentioned above becomes essential. Therefore, the infrastructure design approaches adopted so far get obsolete, which makes paramount the conception and development of advanced circuit modelling techniques. In [11], a thorough analysis of the literature on modelling methods applicable to catenary DC traction systems was conducted. It was discovered from it that nodal analysis-based approaches, mostly conceived via text-based coding, are capable of modelling catenary grids of complicated topology and they also enable flexibility in regulating the number of in-transit vehicles [12], [13], [14], [15], [16], [17]. To summarize, all of the electrical components of the infrastructure are modelled as grid nodes (linked by branches), where voltages and currents are determined iteratively. These approaches allow for quite detailed simulations without requiring excessive time effort, but they lack a userfriendly man-machine interface. In [10], a modular model was proposed in the Matlab-Simulink environment. Due to its modularity, it was still possible to simulate multiple transit vehicles running along contact lines of complicated topology, VOLUME 11, 2023 35023 Authorized licensed use limited to the terms of the applicable license agreement with IEEE. Restrictions apply.  overcoming the limitations of earlier graphical programmingbased modelling techniques for representing the variation in the OCL electrical resistance over time caused by trolleybus motion [18], [19], [20], [21], [22]. As explained in the next section, the methodology introduced in [10] adopts the graphical interface in a way that lends a notable degree of flexibility to the corresponding circuit model. The only issue detected in the aforementioned model concerns an unavoidable spatial discretization of the vehicle's position profiles, imposed by the way the model itself is constructed, the refinement of which may involve considerable computational effort. Given the significance of an intuitive man-machine interface specifically for challenges of this type, this work provides a strategy for improving the precision of modular multi-vehicle motion-based models while minimizing the time required, hence avoiding the use of nodal analysis-based approaches. Bologna's trolleybus network was chosen as a case study, with FSs with simple and, primarily, intricate topologies, with the latter underlining the capabilities of the developed simulation tool.
The discussion is arranged as follows. Section II provides a brief overview on the discrete method (DM) devised in [10]. This is retained useful for clarifying the steps to develop a new continuous modelling method (CM), reported in Section III. Focusing on specific routes within two of the FSs Bologna's traction system is divided into, Section IV includes a comparison under equal conditions between the DM and CM in terms of OCL voltages and currents evaluation, with a view to verify the CM itself. Additionally, by simulating a detail of the network morphology, the advances in model precision attained with the CM are shown. Section V proves the CM's soundness for handling more realistic scenarios. Lastly, the conclusions are drawn in Section VI.

II. USEFUL INSIGHTS ON DISCRETE MULTI-VEHICLE MOTION-BASED CATENARY MODELLING
Similarly to [10], the objective of the model proposed in this paper is the simulation of a trolleybus network's behaviour as a function of the variation in vehicle position over time, also allowing flexibility in selecting the number of in-transit vehicles. For this purpose, detailed analysis of electrical substations is deemed negligible. Indeed, in the literature each TPSS is represented with the Thevenin equivalent circuit [10], [12], [13], [14], [15], [16], [17], [18], [19], [20], [21], [22]. Moreover, the model is decoupled from the representation of the vehicle dynamics. The catenary behaviour is not remarkably affected if the circuit element modelling the vehicle is simply a time-dependent current generator, whose current profile is set as input to the model and may be externally derived through a preceding simulation dedicated just to the vehicle dynamics. This representation of the trolleybus causes the circuit to be linear, which makes it possible to develop the method devised in this work.
As mentioned in the previous section, the issue arising from the representation of the electrical resistance variability in a user-friendly block-based simulator pushed towards a modular modelling approach. Briefly, the absence of a module directly emulating the variable resistor led authors in the field to publish on different expedients that, however, find applications in relatively simple simulation cases (nonmeshed grid topologies), are not suitable for dealing with multi-vehicle scenarios, and sometimes encounter technical criticalities (unphysical conditions or singularities in the circuit analysis [11]). Hence, the model presented in [10] is founded on the idea of keeping the resistive grid configuration invariable while logically allowing for the displacement of the vehicle power absorptions. To this end, the OCL was decomposed into a sufficiently large number of identical LEGO-like blocks (LBs), so that the corresponding contact line resistances are small enough to be reasonably considered fixed. Matlab-Simulink provides a technique to reference a subsystem with the model of a basic LB and create multiple instances of such a subsystem, all of which are synchronized to it. Fig. 2 shows an example of an OCL section. Any OCL section consists of more LBs connected in cascade in such a way that the data to be evaluated propagate in series, hence enabling the spatial distribution of electrical parameters (voltages and currents). At the same time, any type of connection may be considered for the positive and negative wires, which makes possible the model extension to the desired catenary length and morphology.    . 3 shows the circuit model of an LB. The circuit represents a span of the two-wire catenary, with a unique time-dependent current source (symbolised as a controlled source in this article, even though it is independent from external voltages or currents) connected in its midpoint (splitting the wire resistances in four half-span elements R OCLhs ), between the positive and negative poles, and simulating the vehicle. The tricky part was the simulation of the vehicle motion. Essentially, by adopting a suitably programmed logic, only the current generator within the LB whose position corresponds to the position of the vehicle in the real system is switched on. In this way, the equivalent line resistance seen by the turned-on current source changes accordingly, thus imitating the behaviour of an existing scenario. A function henceforth referred to as the weight function is introduced to define whether the trolleybus is travelling within a certain LB. For deciding the activation of the corresponding current generator in the DM, the weight function matrix U is defined as where I gen is a vector including the currents to be assigned to the generator within each LB (N b blocks in total), while I T is a vector containing the currents of the entire trolleybus fleet (N v vehicles in total). Equation (1) agrees with the conditions (2). If the j-th trolleybus (T j , with j = 1, . . . , N v ) is within the b-th LB (LB b , with b = 1, . . . , N b ), then the element u b,j equals one, and the current drawn by T j (i T ,j ) is assigned to i gen,b ; otherwise, u b,j equals zero. The methodology the DM is based on is effective for simulating any network topology and lends itself well to the representation of a multi-vehicle scenario, as the respective electrical absorptions automatically superimpose within (1). However, the precision is closely related to the number of LBs discretizing the OCL section under study. Due to the way the basic block is constructed and in line with (2), it is not possible to enter the actual location of the vehicle inside the span circuit, which is instead rounded and referred to the position of the midpoint current generator. Therefore, the length of the OCL span circuit within each LB determines the discretization step. Obviously, this step could be reduced at will, but with the drawback of an intolerable computational effort. Furthermore, voltage and current measurements along the catenary can only be performed at the LB level (the readings are concatenated block after block to obtain the distribution of the electrical parameters along the entire network), as can be understood by observing the voltmeter and ammeters (available in the library of the utilized software) in Fig. 3, which means that the number of measurements acquired is equal to the number of LBs.
Without abandoning the intuitive graphical user interface and maintaining the model modularity property, a novel approach must be adopted to achieve high-precision without excessive computational effort. Section III illustrates the new methodology to address the above remarks.

III. DEVELOPMENT OF THE PROPOSED CONTINUOUS CATENARY MODEL
The idea of the novel catenary model is still to translate the contact line resistance variability into a current source displacement, but in a continuous way, hence without imposing an a priori discretization strictly related to the number of LBs (that is, to their length). Thanks to the two-port network (2PN) theory, the authors in this paper came up with a methodology to arbitrarily augment the model precision locally, i.e., at any point of interest along the OCL. In the case of a connection of an SESS or an EV supply equipment (EVSE) at a certain DC grid point, it is often not necessary to monitor the voltages and currents of the entire network, as was done in [23] using the DM, but it is sufficient to know these electrical parameters in the neighbouring region, albeit with reasonable precision. The following paragraphs will show the potential of the CM in this respect.

A. UPGRADE FROM STEP-LIKE TO RAMP-LIKE WEIGHT FUNCTIONS
Let us make reference to Fig. 4 for the following considerations. For explanatory purposes, the discussion below is made by assuming only one trolleybus T j travelling along the examined k-th OCL span with single-wire resistance R OCLs,k . For better visualization, the span is as long as the basic block used in the DM (see Fig. 3), so that R OCLs,k = 2R OCLhs . As explained later, the span length can change from block to block in the CM, justifying why the block subscript k (k = 1, . . . , N ′ b , where N ′ b is the number of blocks employed in the CM) is also attributed to the resistance. Furthermore, note that Fig. 4 only depicts the resistance of the positive wire. The reason is that the resistances of the positive and negative poles are assumed to be equal and traversed by the same current, so that the circuit may be simplified by virtually bringing VOLUME 11, 2023 the negative-wire resistance to the positive one (disregarding the negative wire), thus doubling the positive-wire resistive contribution. This hypothesis is not valid in reality because of the presence of the ground current, as well as of a nonuniform distribution of the EB system, or any other electrical connection that generates asymmetry between the positive and negative poles. However, later in this paper, it will be shown that the same assumption can be applied just for the development of the CM, by virtue of the 2PN representation, while the actual electric current discrepancy between positive and negative wires is regulated a posteriori.
Analysing the equivalent electric circuit representing the real system, the in-motion vehicle can be seen as a current generator acting as a slider between the positive and negative contact line wires. As the electrical node moves with the displacement of the source, the corresponding equivalent resistance seen from that point changes with time accordingly. The positive wire resistance is split into two variable contributions by the moving current generator, so that the sum of the respective fractions of the total span resistance always has a unit value. As known, the fractions of resistance are closely related to the length of the corresponding circuit branches, and they are here associated to the following functions: where x k,j is the actual position of T j with respect to the left end of the span under analysis (k-th), l s,k is length of the span itself, i.e., the distance between the neighbouring generators i gen,b−1 and i gen,b (midpoint sources of LB b−1 and LB b , respectively), f k,j and h k,j are complementary functions (their sum is equal to one) relative to the k-th span, considered as further weight functions related to the portions of span resistance generated by T j .
In the DM, the lack of intermediate values of the weight function u b,j in the interval [0,1] forces a discretization dependent on the length of the resistive LB, so that the simulated trolleybus actually ''jumps'' from one block to another, i.e., from the span midpoint to the next one. The reason is that the midpoint current source in Fig. 3 is either turned on (the full vehicle current is assigned to it) or off (see (1), and (2)). In the example of Fig. 4, the vehicle T j belongs to LB b , thus causing u b−1,j to be set to zero and u b,j to gain a unit value, in agreement with (2).
The bottom circuit diagram illustrated in Fig. 4 is the key part of this work. Since, for the reasons already explained, dealing with electrical resistance variability is not handy, the idea of the CM is to apply the same degree of variations (f k,j and h k,j ) of the two resistive contributions characterizing the span to two current generators (left-end source i gen l,k and right-end source i gen r,k , both part of the same k-th basic block) placed at its extrema. Basically, the weight functions follow a linear behaviour, instead of step-like, along the block length. The electrical equivalence of the novel circuit configuration devised here (bottom circuit diagram of Fig. 4) with the real system circuit diagram is proved in Section III-C, thanks to the concept of 2PN.

B. NOVEL BASIC BLOCK
Due to the way the weight functions have been defined for the CM, the configuration of the basic block designed for the DM (refer to Fig. 3) loses its validity, as the continuity of the weight functions themselves would be interrupted if electrical nodes were added between adjacent blocks. Electrical nodes may be generated by the connection of voltage equalizers or feeding points (feeders, reinforcing feeders (RFs), EVSE connection points, BESS connection points, and so forth), which cause the flow of additional currents, thus compromising the subdivision of current contributions between consecutive generators. For this reason, the new basic block has the appearance of the circuit diagram shown in Fig. 5, with two current generators placed laterally representing all the electrical absorption contributions of the vehicles crossing the span under analysis (a proper logic was employed to avoid the simultaneous activation of in-parallel generators at the boundary between adjacent basic blocks).
The DM was not faced with this problem, given that the u function has the form of an indicator function and that all the trolleybuses inside the OCL span were reported to the position of the midpoint current source via a rounding operation (see (2)).
The following paragraphs illustrate the main features of the CM, which lead to its effectiveness. The way the continuous method is structured, it is sufficient to refer to a single generic k-th basic block for the analysis below: from here on, subscript k is omitted.

1) MANAGEMENT OF THE CURRENT DISCREPANCY BETWEEN POSITIVE AND NEGATIVE OCL WIRES
The evaluation of voltages and currents along the stretch of catenary within a basic block can be performed by assuming symmetry between positive and negative wires, by virtue of the 2PN properties. The grounding system, along with any uneven electrical connection with the two wires, justifies the difference between the respective currents.
Let us assume one has already calculated the currents at the desired i-th (i = 1, . . . , M , where M is the number of measurement points) position within each basic block under the condition of wire symmetry. The actual currents flowing along the positive and negative conductors (i OCL+,i and i OCL−,i , respectively) can be expressed as follows: where i com,i can be defined as the common-mode current, i.e., the current already calculated and equal for both wire, while i diff is the differential-mode current, i.e., the contribution to be considered for matching the currents in the real system. Subtracting the equations in (4), one obtains: Making reference to Fig. 5, for any basic block, i OCL+,i and i OCL−,i in (5) assume the current values measured by the corresponding ammeters. Since there are no electrical connections within each block other than the nodes due to the presence of vehicles across the positive and negative contact lines (the vehicles affect both pole currents equally), it is therefore guaranteed that the difference between the values collected by the two ammeters is the same as between the two wires at any point within the base block itself (that is why the subscript i in i diff is omitted).

2) EVALUATION OF VOLTAGES AND CURRENTS ALONG THE OCL
The DM was characterized by the evaluation of catenary voltages and currents at the LB level, i.e., for each LB only a pair of these two quantities was collected via the corresponding voltmeter and ammeter, respectively (actually three quantities, if considering also the ammeter for measuring the negative-pole current, as depicted in Fig. 3).
With the CM, on the other hand, measurements can be taken at any desired number of points within each basic block, allowing arbitrary precision. The number of LBs is thus no longer constrained to the amount of electrical parameters acquisitions, since the latter is now set as input to the model, as explained further on. This allows for a reduction of the amount of LBs to be introduced wherever there are no electrical nodes in the grid. Basically, between any pair of adjacent electrical nodes it is enough to place just one LB of a suitable length, and voltages and currents are calculated inside it at any location decided in the input. Therefore, as the blocks having the configuration reported in Fig. 5 may be endowed with adjustable length, henceforth they are referred to as tunable LBs (TLBs).
Owing to the linearity of the network, the superposition principle was exploited for calculating OCL voltages and currents within each TLB, by summing the responses obtained from each trolleybus considered separately. Being the TLB a 2PN, the electrical parameters determined in each block depend only on the vehicles within the block itself. The contribution to the TLB of each vehicle outside it is concentrated in the across-port voltages measured by the left and right voltmeters, shown in Fig. 5. If no trolleybuses are present inside the TLB, then the common-mode current is flowing only by virtue of the discrepancy between v ms l and v ms r . This current is defined as circulating current: Basics of circuit theory (in particular, the Ohm's law and the current divider rule) are required to obtain the commonmode current and the voltage along the catenary. While in the DM the measurement points were established by the number of LBs, here voltages and currents are calculated at the positions of the vehicles (so that a unique i-th index can be utilized, as one can see in the following paragraphs), which are collected in the position vector injected as input to the model. This is due to the dependence of the electrical parameters on the weight function matrices F and H . Therefore, to add further measurement points, the idea is to enlarge the vehicle position vector in the input (and the corresponding vehicle current vector I T ) with virtual trolleybuses (N vir in total) absorbing zero current. Given N v real vehicles, the total VOLUME 11, 2023 number of measurement points is thus M = N v + N vir . Equations (7) and (8), as shown at the bottom of the page, gather the formulas in matrix form implemented in each TLB for determining respectively OCL common-mode currents (I com ) and voltages (V OCL ). tri l and tri u are functions generating respectively upper triangular and the lower triangular matrices, while tri s u produces a matrix retaining the elements of the original matrix above the first superdiagonal; 1 is an all-ones column vector with M elements. This matrix approach requires the elements to be in a certain order, which would inevitably be compromised owing to the movement of the real vehicles according to their position profile set as input. Hence, the CM requires the matrix reordering at each simulation step. For each TLB, the non-null elements of F and H are those corresponding to the actual or virtual vehicles inside the TLB itself. Hence, there is no need of measurements concatenation as for the DM (refer to Section II), but the contributions of all the TLBs in terms of OCL voltages and currents can be added together without overlapping. If I com,k and V OCL,k are the vectors containing the electrical parameters evaluated in the desired points along the entire catenary only with reference to the k-th TLB, the complete vectors I com and V OCL are obtained as follows: where N ′ b is the number of TLBs constituting the OCL section under analysis. Now, first by adding, then by subtracting the vector with the differential-mode currents (I diff ) determined for all the TLBs, from the I com in (9) the entire distributions respectively of the positive-and negative-wire OCL currents are derived:

3) COMPUTATIONAL EFFORT REDUCTION
The independence of the model precision from the basic block number gives the CM an advantage over the DM in terms of computational effort. However, further steps make such a benefit even sharper. As said already, the matrix calculations summarized in the equations (7), and (8) necessitate the systematic rearrangement of the matrix elements, which affects the simulation time. In addition, we distinguish two cases that explain how in the matrix approach often unnecessary or trivial calculations are made, which can be avoided. Consider the model of a section of a meshed trolleybus network (such as that of the city of Bologna), thus characterised by several TLBs interposed between electrical nodes. In the case in which the measurement points of the electrical parameters are distributed along the entire model of the grid section, in calculating voltages and currents within a given TLB, the functions f i and h i will assume a null value for all the measurement points outside that TLB: the matrices in (7),and (8) therefore tend to be sparse. If, on the other hand, measurement points are concentrated locally in a TLB, all matrix calculations involving the other TLBs are useless as they would not contribute to providing the final voltage and current values expressed in (9), and (10). The methodology adopted in this work thus involves the decomposition of the matrix calculation through the use of nested iterative control structures (for loops) guided by conditional instructions (if statements). The flow chart shown in Fig. 6 helps visualise this. For each TLB, OCL voltages and currents are evaluated only if the measurement point is within the TLB itself, i.e., when f i + h i = 1. If so, the loop control 35028 VOLUME 11, 2023 Authorized licensed use limited to the terms of the applicable license agreement with IEEE. Restrictions apply. variable t is introduced to perform the superposition principle for computing the catenary electrical parameters.

C. CIRCUIT THEORY BEHIND THE MODEL GENERATION
To prove the validity of the translation of contact line resistance variability (real system) into current source variability (proposed model), i.e., to verify if the weight functions reported in (4) can be ascribed to the two side current sources in the TLB, fundamentals of circuit theory are resorted to in this paper. Section III-B1 showed that the discrepancy between the positive and negative wire current can be accounted for even after voltage and current calculations within the basic block. Hence, for demonstrating the electrical equivalence between the OCL span electrical representation and the basic block for the CM, the resistance of the negative wire can be omitted, as illustrated in Fig. 7. For the sake of simplicity, only one vehicle (real or virtual) T i is travelling along the span under study. As already mentioned, any stretch of catenary can be represented as a 2PN, with intoport currents i 1 and i 2 , and across-port voltages v 1 and v 2 . if the two circuit diagrams in Fig. 7 have the same two-port equivalent Norton representation, depicted in Fig. 7, then they are equivalent to each other. The main steps to verify such an electrical equivalence are given below. As known, the conductances of the Norton equivalent circuit (see Fig. 7) are determined by deactivating all the independent sources, i.e., by disconnecting the current generators modelling the trolleybuses. Simply by observing Fig. 7 and remembering that f i + h i = 1, the input and output conductances are derived as g 1,1 = g 2,2 = 1 2R OCLs , equal because of the network symmetry. Analogously, the direct and reverse transfer conductances are expressed as g 1,2 = g 2,1 = − 1 2R OCLs , coincident due to the network reciprocity. Now, the short-circuit currents of the 2PN are calculated in (11), and they are represented by the equivalent independent current generators shown in the circuit diagram of Fig. 7. The current divider rule is adopted for transforming the real system schematic drawing (left side of Fig. 7) into the Norton equivalent circuit, remembering that f i + h i = 1: The considerations above confirm the reasonableness of the basic block employed in the CM, also reported in Fig. 5. Given the linearity of the circuit under study, the analysis can be extended for a multi-vehicle case thanks to the superposition principle, which leads to the generalization of the shortcircuit currents in (12): where M s is the number of vehicles (real plus virtual) in the considered OCL span, corresponding to the number of desired measurement points.

IV. MODEL VERIFICATION AND PRECISION ANALYSIS
In the present section, the CM is verified by comparing it with the DM in terms of voltage and current distribution along two FSs (considered in isolation for the sake of simplicity) that constitute Bologna's trolleybus network. Fig. 8a and 8b show the geographical maps of FS ''Marconi-Carducci'' (FS 1) and FS ''Marconi-Trento e Trieste'' (FS 2), respectively. Both FSs analysed here are bilaterally supplied and characterised by a complex structure (which makes the effectiveness of the comparison more evident), owing to the co-presence of double-bifilar lines and single-bifilar loops, as well as to the inclusion of several RFs scattered along the catenary. In addition, a graphical study of the voltage distribution recorded in a zoomed stretch of FS 1 further highlights the recommended model's gain in precision.

A. VERIFICATION OF THE NOVEL CATENARY MODEL
For the DM-CM comparison to be reasonable, the starting conditions of the simulation must be analogous. First, the total bifilar length was rounded to the nearest multiple of 20 m, in agreement with the choice of a span length of 20 m inside each LB. Thus, the catenary lengths of FS 1 (4541 m) and FS 2 (4372 m) were approximated from to 4540 m and 4380 m, respectively. As a consequence, the position of any electrical element connected to the OCL was rounded to the nearest multiple of 20 m. Secondly, the virtual vehicles representing the measurement points in the CM were placed at the positions of the LBs' midpoint generators. Furthermore, in order to obtain the exact correspondence of the catenary voltages and currents with the two methods, an instant of time was captured in which the real moving trolleybuses were also located at the positions of the LBs' midpoint current sources. The simulation scenarios considered make reference to the traction network data reported in [10], where FS 1 was widely  studied. FS 1 was also examined in [23], where considerations were made regarding its possible modernisation through smart city-oriented technologies. FS 2 was selected as a case study for several system-level analysis [24], [25], [26]. Fig. 9 and 10 depict the voltage and current along the contact line of FS 1 and FS 2 respectively, evaluated both with the DM and the CM at an instant of time during which nine real trolleybuses are positioned at specific points dictated by the aforementioned conditions. For the purpose of this article, it was decided to focus only on the current through the positive wire (i OCL+ ). In addition, the open-circuit voltage at the TPSS terminals was tuned to 800 V (corresponding to the voltage on the AC side of the 12-pulse rectifier, i.e., 590 V). It is evident that the results of the DM and CM coincide perfectly without error, as it should be. This proves that the CM works, i.e., the 2PN concept finds valid application in catenary-powered electric traction grids.

B. GRID MORPHOLOGY EXAMINATION FOR MODEL PRECISION TESTING
It was believed useful to explore the complexity of the morphology of Bologna's trolleybus network to show the capability of the suggested circuit modelling technique in terms of precision. Similar research using the DM was presentedÂ in [10]. The significance of measuring the voltage that the vehicle is exposed to over time as it moves along the FS 1 was recognized for this goal. The graph in Fig. 11 illustrates the catenary voltage trend that one trolleybus (the only one in the FS 1 for this study) sensed as a function of its linearly time-varying position during its route inside the FS. It is more noticeable that the voltage determined through the CM follows the behaviour of the one estimated with the DM if the locations of any topological element (EB, feeder, and RF) are rounded to multiples of 20 m in the CM (refer to the right graph in Fig. 11). In comparison to the 20-m spatial discretization that distinguishes the DM, the continuity of the OCL voltage profile obtained with the CM is highlighted, evident also because of the visible angular points generated by EBs and feeding cables.

V. TOWARDS THE SIMULATION OF REAL SCENARIOS
The results obtained in Section IV relate to deliberately simplified scenarios for the purpose of comparison between the DM and CM. However, to show the applicability of the CM to any realistic operation of a trolleybus system, the catenary model should be adjusted accordingly. The baseline scenario is the one discussed in [26], whose main features are reported below.

A. THE BASELINE SCENARIO FOR ACCURATE SIMULATION
Once the model contains the detailed topology of the DC grid and it is able to simulate the displacement of the vehicles over time along the catenary, it is paramount to characterise the vehicle's electrical consumption to make the simulation itself more truthful. Since urban transport time schedules generally repeat daily, it is reasonable to simulate a one-day operation of the trolleybus network. In [26], the authors devise an approach to refine the vehicle's current and position profiles according to voltage and current measurement data at the output of a TPSS, collected over the course of one day. The vehicle's current draw was customized to its specific route in the analysed FS. Multiple factors affect the vehicles' current profile: not just the time schedule and the map of related bus stops, but also phenomena occurring with different orders of frequency. Only what were referred to as low-frequency phenomena (mainly vehicle acceleration and braking) were modelled, while high-frequency random phenomena (driver's behaviour, pedestrians crossing the road, traffic lights, accidents, and so forth) were neglected owing to their minor impact on the network energy content. An important role is played by the vehicle weight variation throughout the day, which is a very low-frequency contribution. Finally, with the help of a model of a typical trolleybus traction system, the correlation between speed and current consumption of the vehicle versus time in any stretch between two bus stops was possible.
To verify the CM, trolleybuses were treated as current sources/loads in Section IV. However, it is more realistic to consider them as power loads (or sources during the braking phase of conventional trolleybuses), as done in [26] with the DM. Section V-B elaborates on this aspect.

B. MODELLING TROLLEYBUSES AS POWER LOADS
Concerning the DM, as said, voltage and current measurements at a given position are acquired respectively by the voltmeter and ammeter within the relative LB, hence the number of measurement points is strictly related to the number of LBs. Therefore, the logic within each LB must be updated for allowing the current-into-power load modification. Once the power of the vehicles (obtained through the trolleybus model discussed in [26]) is injected as input to the cascaded LBs, within each LB it is divided by the voltage v ms,b measured across the midpoint current source (see Fig 3) in order to get the current i T ,j to be drawn by the midpoint generator itself (refer to (1)). On the other hand, the TLB employed in the CM (see Fig. 5) remains unaltered by the current-into-power load transformation. In such a way, each TLB sees the vehicles inside it as current loads, which keeps the circuit linear and thus allows the application of the superposition principle on which the CM is based, without the need of resorting to the classical iterative methods (such as Newton-Raphson or backward/forward sweep method for non-linear network resolution [17]) that might affect the overall computational effort. This is due to the measurement points (not related to the TLB number) being established in the model input via the virtual trolleybus expedient. The OCL voltages obtained in the model output are fed back to the input to divide the vehicle power profiles, aiming at producing the current vector I T to be fed into the cascaded TLBs (for solving 7, and 8).
The above modifications imply that both the DM and CM require iterative loops within the model to obtain the desired vehicle currents. In the following it is explained how model convergence is guaranteed for the simulation cases of interest. Section V-C reports an illustrative example. Fig. 12 schematises the operation of a vehicle within the DC grid in terms of voltage and current. At each time instant, the vehicle behaves as a constant-power load, hence its current-voltage (I-V) characteristic has a hyperbolic shape. On the other hand, the DC grid may be simply represented with a Thevenin or Norton I-V trans-characteristic (straight line). The latter curve intersects the current and voltage axes at two fundamental points: the Norton current (I N ), i.e., the short-circuit current, and the Thevenin voltage (V TH ), i.e., the open-circuit voltage. Excluding any special case, the two I-V trans-characteristics intersect at two distinct operating points, one stable (O s ) and the other unstable (O u ). Given a certain variation in the voltage across the trolleybus, the condition determining the stability of O s is expressed by (13): Any operation of the DC grid at the vehicle's terminals with a voltage higher than the voltage at O u (i.e., working area on the right of O u ) will ensure convergence to O s . Consequently, by setting the no-load substation voltage (i.e., 800 V) as the initial condition for the grid voltage, the convergence of the simulation is guaranteed: indeed, as this value corresponds to V TH , i.e., the highest voltage one can get from the power FIGURE 13. Relationship between the data time step (t ds ) and the solver sub-step (t ss ).
supply, from Fig. 12 it is clear that the operation is certainly on the right of O u .
Similar considerations apply when the vehicle operates in regenerative braking mode. In this case, being the vehicle electric power negative as injected into the DC grid, the load curve is reflected with respect to the horizontal axis in Fig. 12. The power injection causes the OCL voltage to increase beyond V TH . By prolonging the source characteristic in the fourth quadrant of the I-V plane, another intersection point is generated, which is stable according to (13).
In what follows, we focus on the convergence test for the CM, since the DM has already been tested in the work that led to the publication of [26].

C. SIMULATION CONVERGENCE TEST
The non-linearity of the load power profile depicted in Fig. 12 leads the simulator to perform iterations to achieve a certain vehicle current such that the derived power (P v ) converges to the reference power (P * v ) provided as input. For the purposes of this work, the convergence is verified when the relative error ε r between P v and P * v satisfies the following condition: For a better understanding of the convergence process, a fixed-step solver was selected. If f ds is the data sampling frequency, the solver time step size is such that the respective number of sub-steps (or sub-cycles) occurring in one second (f ss ) is an integer multiple of f ds . To verify the convergence of the precise catenary model, a stepped power profile is applied as input. Generally, a step signal is used as a test signal because the corresponding response of a system gives the information about how quickly the system itself reacts to a sudden change in the input signal. As this reaction is one of the most challenging one can expect from a system, it can be argued that if the simulation converges with the application of a stepped input waveform, it will also converge with any real power profile. Fig. 13 schematises a stepped reference power profile by highlighting the relationship between the data time step and the solver sub-step. The convergence test setting data and results are summarised in Table 1, assuming f ss = 5f ds . The simulation lasted two time cycles, the first one characterized by a positive vehicle power step, while the second one by a negative step with sign reversal. For the sake of this analysis, the trolleybus was hypothesised at a fixed position where it followed a fictitious stepped motion profile representing a worst-case scenario from the simulation convergence viewpoint. Convergence has been verified within  five sub-cycles for any position along the section of trolleybus grid under study. As a precautionary assumption, the number of sub-cycles within each time step has been set such that f ss = 10f ds . In doing this, the time sub-step size has been guaranteed to be less than the system's fastest time variation.
Finally, CM convergence can be also justified in modelling more realistic scenarios, in particular a one-day operation of FS ''Marconi-Trento e Trieste'', already investigated in [26]. Depending on the time schedule of the vehicles running in FS 2, a maximum of five trolleybuses are in operation at any one time. In this case study, only conventional trolleybuses with regenerative braking, i.e., which inject current into the DC grid (hence the vehicle power profile features sign reversal) while braking, are employed.
Based on the criteria for accurate simulation reported in Section V-A, Fig. 14 shows the electric power of a vehicle within FS 2 corresponding to two of its operating cycles, or traction diagrams (refer to [26]). Setting f ss = 10f ds and collecting the data every 10 sub-steps result in an exact matching between the calculated power P v and the reference input power P * v . Simulation convergence is thus ensured.

VI. CONCLUSION
A high-precision model of trolleybus networks was presented in this work. To verify the rationality of the suggested method, the paper provides a direct comparison with a modular discrete-space modelling technique that was previously published in the literature. In contrast to the preceding approach, the number of electrical parameter acquisitions is unbound by the quantity of LEGO-like basic resistive blocks composing the grid model, hence avoiding both the inevitable spatial discretization and an excessive computational effort in the event that a more spatially accurate simulation is required. The two-port network concept showed that it is valid to translate the variability of contact line resistance as a function of vehicle motion (not handy to be modelled, as affirmed in the literature) into that of electric current draw, enabling model-independent arbitrary precision in determining the distribution of network voltages and currents. The line parameters between the two compared methods match perfectly when analogous conditions are imposed, but when simulating the morphology in detail, the precise method exalts the continuity of the voltage profile as opposed to a stepped behaviour attributable to the discrete method, according to a graphical analysis of a section of Bologna's trolleybus grid.
Finally, the applicability of the continuous approach to any realistic operation of a trolleybus system was shown. To face the network non-linearity generated by assuming vehicles as power sources/loads, there is no need to resort to conventional iterative methods (more computationally expensive), but the stability of the model is automatically verified for the DC grids of our interest, entailing the convergence of the simulation exclusively via the iterations performed by the software solver.