**3. Calculation of GIC**

#### **3.1 Power networks**

34 Space Science

described by a surface impedance (see equation (9)). The contribution to the total geoelectromagnetic disturbance at the Earth's surface from (secondary) currents in the Earth is an upward-propagating reflected wave. This kind of a model is already included in the basic paper of magnetotellurics by Cagniard (1953). It is necessary to emphasise that the frequencies involved in geoelectromagnetic studies are typically in the mHz range and at least below 1 Hz, i.e. so small that the displacement currents (= t can practically always be neglected. It means that "geoelectromagnetic plane waves" are actually not "waves", so the terminology generally used in geoelectromagnetics is not completely correct

*e i* 

) is associated with the time derivative of *Bx*(*t*), equation (10) can be

(11)

0

Equation (10) shows that there is a 45-degree (/4-radian) phase shift between the geoelectric and geomagnetic fields. We also see that an increase of the angular frequency and a decrease of the Earth's conductivity enhance the geoelectric field with respect to the

inverse-Fourier transformed to give the following time domain convolution integral

0 0 <sup>1</sup> ( ) ( ) *<sup>y</sup> gt u E t du u*

where the time derivative of *Bx*(*t*) is denoted by *g*(*t*). The derivation of equation (10) makes use of the neglect of the displacement currents, which thus also affects equation (11). If the displacement currents are included the kernel function convolved with *g*(*t*) in formula (11) is more complicated containing the Bessel function of the zero order (Pirjola, 1982). Equation (11) is in agreement with the causality, i.e. *Ey*(*t*) at the time *t* only depends on earlier values of *g*(*t*). The square root of the lag time *u* in the denominator means that the influence of a value of *g*(*t-u*) on *Ey*(*t*) decreases with increasing *u*. As indicated in Section 2.1, equation (11) shows that the relation of the geoelectric field (and GIC) with the geomagnetic time derivative is not simple, like for example a proportionality, not even in the present plane-

Similarly to the inverse-Fourier transform of equation (10) leading to (11) in the time domain, we can inverse-Fourier transform equation (9) to get a convolution relation between *Ey* and *Bx* in the time domain, or as above in equation (11), perhaps rather between *Ey* and d*Bx*/d*t*.

The engineering part of GIC modelling utilises the horizontal geoelectric field to be provided by the geophysical part and also needs the knowledge of the topology,

) at the Earth's surface is related to the

<sup>4</sup> *Bx* (10)

. It is easy to show that the

and consider a

) by the following

Let us assume now that the Earth is uniform with the conductivity

perpendicular horizontal geomagnetic variation component *Bx* = *Bx*(

*Ey*

harmonic time dependence with the angular frequency

horizontal geoelectric field component *Ey* = *Ey*(

in this respect.

equation (e.g. Pirjola, 1982)

geomagnetic field.

*Bx*(

(Cagniard, 1953; Pirjola, 1982)

wave and uniform-Earth situation.

Noting that i

For investigating the engineering part of GIC modelling in the case of a power system, we consider a network of conductors with *N* discrete nodes, called stations and earthed by the resistances *Re,i* (*i* = 1,…,*N*). Let us assume that the network is impacted by a horizontal geoelectric field **E**, which implies the flow of geomagnetically induced currents (GIC). Lehtinen and Pirjola (1985) derive a formula for the *N* 1 column matrix **Ie** that includes the currents *Ie,m* (*m* = 1,…,*N*) called earthing currents or earthing GIC and flowing between the network and the Earth as follows

$$\mathbf{I}\_e = (\mathbf{1} + \mathbf{Y}\_\mathbf{n} \mathbf{Z}\_e)^{-1} \mathbf{J}\_e \tag{12}$$

The current is defined to be positive when it flows from the network to the Earth and negative when it flows from the Earth to the network. The symbol **1** denotes the *N N* unit identity matrix. The *N N* earthing impedance matrix **Ze** and the *N N* network admittance matrix **Yn** as well as the *N* 1 column matrix **Je** are explained below.

The definition of **Ze** states that multiplying the earthing current matrix **Ie** by **Ze** gives the voltages between the earthing points and a remote Earth that are related to the flow of the currents *Ie,m* (*m* = 1,…,*N*). Thus, expressing the voltages by an *N* 1 column matrix **U**, we have

$$\mathbf{U} = \mathbf{Z}\_{\mathbf{e}} \mathbf{I}\_{\mathbf{e}} \tag{13}$$

Utilising the reciprocity theorem, **Ze** can be shown to be a symmetric matrix. The diagonal elements of **Ze** equal the earthing resistances of the stations. If the distances of the stations

Geomagnetically Induced Currents as Ground Effects of Space Weather 37

resistances of possible reactors or any resistors in the earthing leads of transformer neutrals,

As Mäkinen (1993) and Pirjola (2005) present, special modelling techniques are needed for transformers when two different voltage levels are considered in the calculation of GIC in a power network. Moreover, the treatment of autotransformers differs from the case of full two-winding transformers. It is worth mentioning that recent, still unpublished, investigations indicate that the methods and equations described by Mäkinen (1993) and

Fig. 3. Finnish 400-kV electric power transmission grid in its configuration valid in October 1978 to November 1979 (Pirjola & Lehtinen, 1985; Pirjola, 2005; Pirjola, 2009). This network

Pirjola (2009) introduces an old version of Finland's 400 kV power network as a "test model" for GIC computation algorithms and programs. The network is shown in Figure 3. It consists of 17 stations and 19 transmission lines, so that it is complex enough to reveal essential features included in GIC computations but not too large to unnecessarily complicate the calculations and analyses. Pirjola (2009) provides the (Cartesian) north and east coordinates of the locations of the stations included in the test model and numbered from 1 to 17. The (total) earthing resistances of the stations as well as the resistances of the transmission lines are also given. Regarding the use of the network as a test model, Pirjola (2009) presents GIC values at all stations and in every transmission line that a uniform (geographically) eastward or a uniform (geographically) northward geoelectric field of 1

can be used as a "test model" for GIC calculation algorithms and programs.

with all these resistances connected in series.

Pirjola (2005) are unnecessarily complicated, though correct.

are large enough, the influence of the current *Ie,m* at one station on the voltages at other stations is negligible, and then the off-diagonal elements of **Ze** are zero (see Pirjola, 2008).

The matrix **Yn** is defined by

$$(\text{i} \neq m): Y\_{n,im} = -\frac{1}{R\_{n,im}}, \ (\text{i} = m): Y\_{n,im} = \sum\_{k=1, k \neq i}^{N} \frac{1}{R\_{n,ik}} \tag{14}$$

where *Rn,im* is the resistance of the conductor between stations *i* and *m* (*i*, *m* = 1,…,*N*). (If stations *i* and *m* are not directly connected by a conductor, *Rn,im* naturally gets an infinite value.) It is seen from equation (14) that **Yn** is symmetric.

The elements *Je,m* (*m* = 1,…,*N*) of the column matrix **Je** are given by

$$J\_{e,m} = \sum\_{i=1, i \neq m}^{N} \frac{V\_{im}}{R\_{n,im}} \tag{15}$$

The geovoltage *Vim* is produced by the horizontal geoelectric field **E** along the path defined by the conductor line from station *i* to station *m* (*i*, *m* = 1,…,*N*), i.e.

$$V\_m = \bigcap\_i^m \mathbf{E} \cdot d\mathbf{s} \tag{16}$$

Generally, the horizontal geoelectric field is rotational, i.e. the vertical component of the curl of **E** is not zero. As seen from equation (3), this is the case when the time derivative of the vertical component of the magnetic field differs from zero. Consequently, the integral in equation (16) is path-dependent. Thus, as indicated, the integration route must follow the conductor between *i* and *m* (Boteler and Pirjola, 1998b; Pirjola, 2000). Equations (15) and (16) show that the column matrix **Je** involves the contribution from the geoelectric field **E** to equation (12). Note that **Ie** and **Je** are equal in the case of perfect earthings, i.e. when **Ze** = 0.

Pirjola (2007) presents an alternative, but equivalent, expression for the matrix **Ie**, which includes the geovoltages more explicitly. It makes use of the total admittance matrix that is defined to be the sum of **Yn** and the inverse of **Ze**. However, the application of this alternative expression obviously does not produce any noticeable advantages, for example, in numerical computations.

So far, we have discussed GIC flowing between the network and the Earth. There is also a simple expression for GIC in the conductors between the nodes (e.g. Pirjola, 2007). However, we ignore it and concentrate only on earthing GIC because, in the case of a power system, GIC to (or from) the Earth constitute the harmful currents in transformer windings and they are usually measured (see Section 2.1).

When GIC in a power network are calculated, the three phases are usually treated as one circuit element. The resistance of the element is then one third of that of a single phase, and the GIC flowing in the element is three times the current in a single phase. Moreover, the earthing resistances (which might be called the *total* earthing resistances) are convenient to be assumed to include the actual earthing resistances, the transformer resistances and the 36 Space Science

are large enough, the influence of the current *Ie,m* at one station on the voltages at other stations is negligible, and then the off-diagonal elements of **Ze** are zero (see Pirjola, 2008).

, ,

*i mY i mY*

,

value.) It is seen from equation (14) that **Yn** is symmetric.

The elements *Je,m* (*m* = 1,…,*N*) of the column matrix **Je** are given by

by the conductor line from station *i* to station *m* (*i*, *m* = 1,…,*N*), i.e.

*n im n im*

1 1 ( ): , ( ):

where *Rn,im* is the resistance of the conductor between stations *i* and *m* (*i*, *m* = 1,…,*N*). (If stations *i* and *m* are not directly connected by a conductor, *Rn,im* naturally gets an infinite

> *N im e m i im n im*

The geovoltage *Vim* is produced by the horizontal geoelectric field **E** along the path defined

*m*

Generally, the horizontal geoelectric field is rotational, i.e. the vertical component of the curl of **E** is not zero. As seen from equation (3), this is the case when the time derivative of the vertical component of the magnetic field differs from zero. Consequently, the integral in equation (16) is path-dependent. Thus, as indicated, the integration route must follow the conductor between *i* and *m* (Boteler and Pirjola, 1998b; Pirjola, 2000). Equations (15) and (16) show that the column matrix **Je** involves the contribution from the geoelectric field **E** to equation (12). Note that **Ie** and **Je** are equal in the case of perfect earthings, i.e.

Pirjola (2007) presents an alternative, but equivalent, expression for the matrix **Ie**, which includes the geovoltages more explicitly. It makes use of the total admittance matrix that is defined to be the sum of **Yn** and the inverse of **Ze**. However, the application of this alternative expression obviously does not produce any noticeable advantages, for example,

So far, we have discussed GIC flowing between the network and the Earth. There is also a simple expression for GIC in the conductors between the nodes (e.g. Pirjola, 2007). However, we ignore it and concentrate only on earthing GIC because, in the case of a power system, GIC to (or from) the Earth constitute the harmful currents in transformer windings and they

When GIC in a power network are calculated, the three phases are usually treated as one circuit element. The resistance of the element is then one third of that of a single phase, and the GIC flowing in the element is three times the current in a single phase. Moreover, the earthing resistances (which might be called the *total* earthing resistances) are convenient to be assumed to include the actual earthing resistances, the transformer resistances and the

*m i*

*<sup>V</sup> <sup>J</sup>*

1, ,

, , 1,

*R R*

*n im k ki n ik*

(14)

*N*

*<sup>R</sup>* (15)

*V d* **E s** (16)

The matrix **Yn** is defined by

when **Ze** = 0.

in numerical computations.

are usually measured (see Section 2.1).

resistances of possible reactors or any resistors in the earthing leads of transformer neutrals, with all these resistances connected in series.

As Mäkinen (1993) and Pirjola (2005) present, special modelling techniques are needed for transformers when two different voltage levels are considered in the calculation of GIC in a power network. Moreover, the treatment of autotransformers differs from the case of full two-winding transformers. It is worth mentioning that recent, still unpublished, investigations indicate that the methods and equations described by Mäkinen (1993) and Pirjola (2005) are unnecessarily complicated, though correct.

Fig. 3. Finnish 400-kV electric power transmission grid in its configuration valid in October 1978 to November 1979 (Pirjola & Lehtinen, 1985; Pirjola, 2005; Pirjola, 2009). This network can be used as a "test model" for GIC calculation algorithms and programs.

Pirjola (2009) introduces an old version of Finland's 400 kV power network as a "test model" for GIC computation algorithms and programs. The network is shown in Figure 3. It consists of 17 stations and 19 transmission lines, so that it is complex enough to reveal essential features included in GIC computations but not too large to unnecessarily complicate the calculations and analyses. Pirjola (2009) provides the (Cartesian) north and east coordinates of the locations of the stations included in the test model and numbered from 1 to 17. The (total) earthing resistances of the stations as well as the resistances of the transmission lines are also given. Regarding the use of the network as a test model, Pirjola (2009) presents GIC values at all stations and in every transmission line that a uniform (geographically) eastward or a uniform (geographically) northward geoelectric field of 1

Geomagnetically Induced Currents as Ground Effects of Space Weather 39

exact values depending on the radius of the pipeline at the particular section. Thus, the

Concerning GIC and pipe-to-soil voltages, electrically long pipelines (*length >> adjustment distance*) behave differently than electrically short pipes (*length << adjustment distance*). For a long pipeline, the voltage decays exponentially at a distance comparable to the adjustment distance, when moving from either end of the pipeline towards the centre where it is practically zero. In the same areas near the ends, GIC flows between the pipe and the soil, and in the central parts the current along the pipe is spatially constant. For a short pipeline, the voltage changes linearly along the pipeline, and the current along the pipe is small when the ends are insulated from the Earth. Consequently, central parts of a long pipeline are not critical regarding corrosion problems due to GIC. On the other hand, the pipe-to-soil voltages remain smaller in short pipelines than the end voltages of a long pipeline. Therefore, it might sometimes be reasonable to install insulating flanges in series in a

The handling of inhomogeneities associated with a pipeline, such as bends, changes of the pipeline material or of the pipeline size, and branches of the pipeline network, can be accomplished by applying Thévenin's theorem. It expresses an equivalent voltage source and an impedance that describe an external circuit at the terminals of a network section considered. Thus at a pipeline inhomogeneity, Thévenin's voltage and impedance of a section terminating at the inhomogeneity have to be calculated when considering GIC and pipe-to-soil voltages in another section. In this way it is possible to go through a whole pipeline network from section to section. Branches require a special treatment since, for a particular section, Thévenin's components of the other sections are connected in parallel

Geomagnetically induced currents (GIC) are ground effects of space weather, which is associated with a complex chain of phenomena extending from processes in the Sun to GIC in technological networks. In general, GIC are a possible source of problems to the network. Thus research of GIC is practically important, but it also has scientific significance because a ground-based network carrying GIC can be regarded as a huge antenna that collects

The first GIC observations date back to early telegraph systems more than 150 years ago. Today power systems constitute the most important target of GIC research. The problems in power networks result from half-cycle saturation of transformers created by dc-like GIC. In the worst cases, large areas may experience a blackout due to GIC and transformers can be permanently damaged. The most significant GIC event so far is the blackout in Québec, Canada, for several hours during a large geomagnetic storm in March 1989. Another well-known GIC event caused a blackout in southern Sweden during the so-called "Halloween storm" at the end of October 2003. In 2011 to 2014, an EU-funded project is going on in which GIC in the entire European high-voltage system are

In this paper, we discuss the techniques readily available for calculating GIC values in power networks and pipelines. Future research efforts should be focussed on the application

adjustment distance has values from about 20 km to about 60 km.

pipeline to break it into shorter sections.

information from processes in space and within the Earth.

(Pulkkinen et al., 2001b).

**4. Conclusion** 

considered.

V/km creates. Here it should be noted that, by using a linear superposition, these data enable the computation of GIC due to any uniform horizontal geoelectric field impacting the network. For some additional details about the test model, Pirjola (2009) is referred to.

The Rauma 400 kV station discussed in Section 2.1 and especially in Figure 1 was not yet included in the Finnish 400 kV system in October 1978 to November 1979. It is located in the line between stations 4 and 7 quite near station 7. Note that geographic and geomagnetic latitudes are far from being parallel in Finland and that the former are higher. In North America, the geographic latitudes are lower than the geomagnetic since the north geomagnetic pole is on the American side of the geographic pole.

#### **3.2 Pipelines**

Pirjola & Lehtinen (1985) present theoretical computations of GIC for the Finnish natural gas pipeline by using the matrix formalism discussed in Section 3.1 and appropriate to a discretely-earthed network, such as a power system. The assumption is included in the treatment by Pirjola & Lehtinen (1985) that the insulating coating of the pipeline is ideal and perfect with a zero conductivity and that the pipeline is earthed at the cathodic protection stations. This approximation should not be considered very good because the large surface area of the pipe makes the pipeline continuously earthed in practice even though the conductivity of the coating material is very small. Furthermore, the CP stations do not constitute normal earthings since the current can only go to the Earth there (to return from the ground to the pipe elsewhere).

Viljanen (1989) presents a GIC study about the Finnish natural gas pipeline based on the simplified assumption that the pipeline is an infinitely long multi-layered cylindrical structure in a homogeneous medium. The model is in agreement with the continuous earthing but it is otherwise too much idealised. According to this model, GIC flowing along the pipeline may reach values of hundreds of amps, which are clearly larger than those measured in Finland (Section 2.1). Although Viljanen (1989) also estimates the effects of a horizontal change of the Earth's conductivity, and of a bend of the pipeline, the treatment is not yet complete for a real pipeline network. A significant improvement in theoretical modelling is shown by Boteler (1997) by incorporating the distributed-source transmission line theory into pipeline-GIC calculations. In fact, the applicability of the DSTL theory is already considered by Boteler and Cookson (1986). An extension is provided by Pulkkinen et al. (2001b) as they also treat branches of a pipeline network.

In the DSTL theory, the pipeline is considered a transmission line containing a series impedance (or resistance due to the dc treatment) *Z* determined by the properties of the pipeline steel and a parallel admittance *Y* associated with the resistivity of the coating. The geoelectric field affecting the pipeline forms the distributed source. An important parameter, called the adjustment distance, is the inverse of the propagation constant defined by

$$
\gamma = \sqrt{ZY} \tag{17}
$$

Typical values of the adjustment distance of a real pipeline are tens of kilometres. For the Finnish natural gas pipeline, *Z* = 5…9·10–3 km–1 and *Y* = 5·10–2…0.25 –1 km–1 with the exact values depending on the radius of the pipeline at the particular section. Thus, the adjustment distance has values from about 20 km to about 60 km.

Concerning GIC and pipe-to-soil voltages, electrically long pipelines (*length >> adjustment distance*) behave differently than electrically short pipes (*length << adjustment distance*). For a long pipeline, the voltage decays exponentially at a distance comparable to the adjustment distance, when moving from either end of the pipeline towards the centre where it is practically zero. In the same areas near the ends, GIC flows between the pipe and the soil, and in the central parts the current along the pipe is spatially constant. For a short pipeline, the voltage changes linearly along the pipeline, and the current along the pipe is small when the ends are insulated from the Earth. Consequently, central parts of a long pipeline are not critical regarding corrosion problems due to GIC. On the other hand, the pipe-to-soil voltages remain smaller in short pipelines than the end voltages of a long pipeline. Therefore, it might sometimes be reasonable to install insulating flanges in series in a pipeline to break it into shorter sections.

The handling of inhomogeneities associated with a pipeline, such as bends, changes of the pipeline material or of the pipeline size, and branches of the pipeline network, can be accomplished by applying Thévenin's theorem. It expresses an equivalent voltage source and an impedance that describe an external circuit at the terminals of a network section considered. Thus at a pipeline inhomogeneity, Thévenin's voltage and impedance of a section terminating at the inhomogeneity have to be calculated when considering GIC and pipe-to-soil voltages in another section. In this way it is possible to go through a whole pipeline network from section to section. Branches require a special treatment since, for a particular section, Thévenin's components of the other sections are connected in parallel (Pulkkinen et al., 2001b).
