The basic challenge in solving our problem is
that it is defined in an unbounded domain; therefore, certain simplifications
are necessary to formulate it for numerical computations.
First, we restrict our radial domain to the region between the inner
radius 
and the outer radius 
. We introduce
a new dimensionless radial coordinate 
where 
is an arbitrary unit of length. At the inner and outer radii we
put 
and 
.
Second, we introduce a new vertical
coordinate, 
for equations (5)-(6) and
for equation (7). This transformation,
adopted from the method originally used by Jeeps (1975) [7] in the context of spherical dynamos,
maps equations (5)-(7) onto a single, finite
rectangular spatial domain 
, 
. The continuity of a magnetic field
on the surfaces of the disc is taken care of by means of the appropriate
boundary conditions at 
. The condition that A vanishes at
(the total potential approaches 
) translates
into 
. We use an equispaced grid in z
that gives values of the magnetic field within a disc at vertical locations
spaced uniformly with respect to the local disc's thickness.
Outside the disc magnetic field is defined at nonuniformly
spaced vertical locations that are concentrated towards the surface of the disc. In an accretion disc, the desired radial spatial resolution varies with
the radius. Thus, an equispaced radial grid is impractical, as it should
be set by the smallest desired resolution. Instead, we introduce the
new radial variable x, such that 
.
We use an equispaced grid in x over its entire domain 
, which gives values of the magnetic field
at nonuniform radial locations that are concentrated towards the star.
The set of equations (5)-(7) is transformed through
the map 
and discretised using
finite differences second-order accurate in space and first-order
accurate in time. The set of finite differences equations resulting from
equations (5)-(6), which describe a magnetic
field inside the disc, is solved by the Euler method
adopting the Courant condition to ensure numerical stability. Once the
inside equations have been advanced one time step, the outside field is
calculated from the discretised form of equation (7) using a
simple relaxation scheme and enforcing the boundary conditions at
the disc's surfaces.