Additional possible analyses

class magic.bLayers.BLayers(iplot=False, quiet=False)[source]

This class allows to determine the viscous and thermal boundary layers using several classical methods (slope method, peak values, dissipation rates, etc.). It uses the following files:

This function can thus only be used when both powerR.TAG and bLayersR.TAG exist in the working directory.

Warning

This function works well as long as rigid boundaries and fixed temperature boundary conditions are employed. Other combination of boundary conditions (fixed fluxes and/or stress-free) might give wrong results, since boundary layers become awkward to define in that case.

Since this function is supposed to use time-averaged quantities, the usual procedure is first to define the initial averaging time using AvgField: (this needs to be done only once)

>>> a = AvgField(tstart=2.58)

Once the tInitAvg file exists, the boundary layer calculation can be done:

>>> bl = BLayers(iplot=True)
>>> # print the formatted output
>>> print(bl)
__init__(iplot=False, quiet=False)[source]
Parameters:
  • iplot (bool) – display the result when set to True (default False)

  • quiet (bool) – less verbose when set to True (default is False)

__str__()[source]

Formatted output

plot()[source]

Plotting function

magic.bLayers.getAccuratePeaks(rad, uh, uhTop, uhBot, ri, ro)[source]

This functions performs a spline extrapolation around the maxima of the input array uh to define a more accurate location of the boundary layer.

Parameters:
  • rad (numpy.ndarray) – radius

  • uh (numpy.ndarray) – the horizontal velocity profile

  • uhTop (float) – first peak value of uh close to the outer boundary

  • uhBot (float) – first peak value of uh close to the inner boundary

  • ri (float) – the inner core radius

  • ro (float) – the outer core radius

Returns:

four floats: thickness of the bottom boundary layer, thickness of the top boundary layer, extrapolated value of uh at the bottom boundary layer, extrapolated value of uh at the top boundary layer

Return type:

list

magic.bLayers.getMaxima(field)[source]

This function determines the local maxima of the input array field

Parameters:

field (numpy.ndarray) – the input array

Returns:

a list containing the indices of the local maxima

Return type:

list

magic.bLayers.integBotTop(rad, field, ri, ro, lambdai, lambdao, normed=False)[source]

This function evaluates the radial integral of the input array field in the bottom and top boundary layers separately.

Parameters:
  • rad (numpy.ndarray) – radius

  • field (numpy.ndarray) – the input radial profile

  • ri (float) – the inner core radius

  • ro (float) – the outer core radius

  • lambdai (float) – thickness of the inner boundary layer

  • lambdao (float) – thickness of the outer boundary layer

  • normed (bool) – when set to True, the outputs are normalised by the volumes of the boundary layers. In that case, the outputs are volume-averaged quantities.

Returns:

two floats that contains the bottom and top boundary layers integrations (integBot, integTop)

Return type:

list

magic.bLayers.integBulkBc(rad, field, ri, ro, lambdai, lambdao, normed=False)[source]

This function evaluates the radial integral of the input array field in the boundary layer and in the bulk separately.

Parameters:
  • rad (numpy.ndarray) – radius

  • field (numpy.ndarray) – the input radial profile

  • ri (float) – the inner core radius

  • ro (float) – the outer core radius

  • lambdai (float) – thickness of the inner boundary layer

  • lambdao (float) – thickness of the outer boundary layer

  • normed (bool) – when set to True, the outputs are normalised by the volumes of the boundary layers and the fluid bulk, respectively. In that case, the outputs are volume-averaged quantities.

Returns:

two floats that contains the boundary layer and the bulk integrations (integBc, integBulk)

Return type:

list

class magic.ThetaHeat(iplot=False, angle=10, pickleName='thHeat.pickle', quiet=False)[source]

This class allows to conduct some analysis of the latitudinal variation of the heat transfer. It relies on the movie files ATmov.TAG and AHF_mov.TAG. As it’s a bit time-consuming, the calculations are stored in a python.pickle file to quicken future usage of the data.

Since this function is supposed to use time-averaged quantities, the usual procedure is first to define the initial averaging time using AvgField: (this needs to be done only once)

>>> a = AvgField(tstart=2.58)

Once the tInitAvg file exists, the latitudinal heat transfer analysis can be done using:

>>> # For chunk-averages over 10^\degree in the polar and equatorial regions.
>>> th = ThetaHeat(angle=10)
>>> # Formatted output
>>> print(th)
__init__(iplot=False, angle=10, pickleName='thHeat.pickle', quiet=False)[source]
Parameters:
  • iplot (bool) – a boolean to toggle the plots on/off

  • angle (float) – the integration angle in degrees

  • quiet (bool) – a boolean to switch on/off verbose outputs

PickleName:

calculations a

__str__()[source]

Formatted outputs

>>> th = ThetaHeat()
>>> print(th)
plot()[source]

Plotting function

class magic.cyl.Cyl(ivar=1, datadir='.', ns=None)[source]

This class allows to extrapolate a given graphic file on a cylindrical grid. Once done, the extrapolated file is stored in a python.pickle file. It is then possible to display 2-D cuts of the extrapolated arrays (radial cuts, phi-averages, equatorial cuts, z-averages and phi-slices)

Warning

This process is actually very demanding and it might take a lot of time to extrapolate the G_#.TAG file. Be careful when choosing the input value of ns!

>>> # Extrapolate the G file to the cylindrical grid (ns=128, nz=2*ns)
>>> c = Cyl(ivar=1, ns=128)
>>> # Radial cut of v_r
>>> c.surf(field='vr', r=0.8)
>>> # Vertical average of B_\phi
>>> c.avgz(field='Bphi', cm='seismic', levels=33)
>>> # Azimuthal average of v_\phi
>>> c.avg(field='Bphi')
>>> # Equatorial cut of of v_theta
>>> c.equat(field='vtheta')
__init__(ivar=1, datadir='.', ns=None)[source]
Parameters:
  • ivar (int) – the number of the Graphic file

  • datadir (str) – working directory

  • ns (int) – number of grid points in the radial direction

avg(field='Bphi', levels=16, cm='RdYlBu_r', normed=True, vmax=None, vmin=None)[source]

Plot the azimutal average of a given field.

>>> c = Cyl(ns=65)
>>> # Azimuthal average of B_r
>>> c.avg(field='Br', cm='seismic', levels=33)
Parameters:
  • field (str) – name of the input field

  • levels (int) – number of contour levels

  • cm (str) – name of the color map

  • normed (bool) – when set to True, the contours are normalised fro -max(field), max(field)

  • vmin (float) – truncate the contour levels to values > vmin

  • vmax (float) – truncate the contour levels to values < vmax

avgz(field='vs', levels=16, cm='RdYlBu_r', normed=True, vmin=None, vmax=None, avg=False)[source]

Plot the vertical average of a given field.

>>> c = Cyl(ns=65)
>>> # Vertical average of v_s
>>> c.avg(field='vs', cm='seismic', levels=33)
Parameters:
  • field (str) – name of the input field

  • levels (int) – number of contour levels

  • cm (str) – name of the color map

  • normed (bool) – when set to True, the contours are normalised fro -max(field), max(field)

  • vmin (float) – truncate the contour levels to values > vmin

  • vmax (float) – truncate the contour levels to values < vmax

  • avg (bool) – when set to True, an additional figure with the phi-average profile is also displayed

equat(field='vs', levels=16, cm='RdYlBu_r', normed=True, vmax=None, vmin=None)[source]

Plot an input field in the equatorial plane.

>>> c = Cyl(ns=65)
>>> # Equatorial cut of v_\phi
>>> c.equat(field='vphi', cm='seismic', levels=33)
Parameters:
  • field (str) – name of the input field

  • levels (int) – number of contour levels

  • cm (str) – name of the color map

  • normed (bool) – when set to True, the contours are normalised fro -max(field), max(field)

  • vmin (float) – truncate the contour levels to values > vmin

  • vmax (float) – truncate the contour levels to values < vmax

slice(field='Bphi', lon_0=0.0, levels=16, cm='RdYlBu_r', normed=True)[source]

Plot an azimuthal slice of a given field.

>>> c = Cyl(ns=65)
>>> # Slices of v_r at 30 and 60 degrees
>>> c.slice(field='vr', lon_0=[30, 60])
Parameters:
  • field (str) – name of the input field

  • lon_0 (float or list) – the longitude of the slice in degrees, or a list of longitudes

  • levels (int) – number of contour levels

  • cm (str) – name of the color map

  • normed (bool) – when set to True, the contours are normalised fro -max(field), max(field)

surf(field='Bphi', r=0.85, vmin=None, vmax=None, levels=16, cm='RdYlBu_r', normed=True, figsize=None)[source]

Plot the surface distribution of an input field at a given input radius (normalised by the outer boundary radius).

>>> c = Cyl(ns=65)
>>> # Surface plot of B_\phi from -10 to 10
>>> c.surf(field='Bphi', r=0.6, vmin=-10, vmax=10, levels=65)
Parameters:
  • field (str) – name of the input field

  • r (float) – radial level (normalised to the outer boundary radius)

  • levels (int) – number of contour levels

  • cm (str) – name of the color map

  • normed (bool) – when set to True, the contours are normalised fro -max(field), max(field)

  • vmin (float) – truncate the contour levels to values > vmin

  • vmax (float) – truncate the contour levels to values < vmax

magic.cyl.sph2cyl(g, ns=None, nz=None)[source]

This function interpolates the three flow (or magnetic field) component of a G_#.TAG file on a cylindrical grid of size (ns, nz).

Warning

This might be really slow!

Parameters:
  • g (magic.MagicGraph) – input graphic output file

  • ns (int) – number of grid points in the radial direction

  • nz (int) – number of grid points in the vertical direction

Returns:

a python tuple of five numpy.ndarray (S,Z,vs,vp_cyl,vz). S[nz,ns] is a meshgrid that contains the radial coordinate. Z[nz,ns] is a meshgrid that contains the vertical coordinate. vs[nz,ns] is the radial component of the velocity (or magnetic field), vp_cyl[nz,ns] the azimuthal component and vz[nz,ns] the vertical component.

Return type:

tuple

magic.cyl.sph2cyl_plane(data, rad, ns)[source]

This function extrapolates a phi-slice of a spherical shell on a cylindrical grid

>>> # Read G_1.test
>>> gr = MagicGraph(ivar=1, tag='test')
>>> # phi-average v_\phi and s
>>> vpm = gr.vphi.mean(axis=0)
>>> sm = gr.entropy.mean(axis=0)
>>> # Interpolate on a cylindrical grid
>>> Z, S, outputs = sph2cyl_plane([vpm, sm], gr.radius, 512, 1024)
>>> vpm_cyl, sm_cyl = outputs
Parameters:
  • data (list(numpy.ndarray)) – a list of 2-D arrays [(ntheta, nr), (ntheta, nr), …]

  • rad (numpy.ndarray) – radius

  • ns (int) – number of grid points in s direction

Returns:

a python tuple that contains two numpy.ndarray and a list (S,Z,output). S[nz,ns] is a meshgrid that contains the radial coordinate. Z[nz,ns] is a meshgrid that contains the vertical coordinate. output=[arr1[nz,ns], …, arrN[nz,ns]] is a list of the interpolated array on the cylindrical grid.

Return type:

tuple

magic.cyl.zavg(input, radius, ns, minc, save=True, filename='vp.pickle', normed=True, colat=None)[source]

This function computes a z-integration of a list of input arrays (on the spherical grid). This works well for 2-D (phi-slice) arrays. In case of 3-D arrays, only one element is allowed (too demanding otherwise).

Parameters:
  • input (list(numpy.ndarray)) – a list of 2-D or 3-D arrays

  • radius (numpy.ndarray) – spherical radius

  • ns (int) – radial resolution of the cylindrical grid (nz=2*ns)

  • minc (int) – azimuthal symmetry

  • save (bool) – a boolean to specify if one wants to save the outputs into a pickle (default is True)

  • filename (str) – name of the output pickle when save=True

  • normed (bool) – a boolean to specify if ones wants to simply integrate over z or compute a z-average (default is True: average)

  • colat (numpy.ndarray) – an optional array containing the colatitudes

Returns:

a python tuple that contains two numpy.ndarray and a list (height,cylRad,output) height[ns] is the height of the spherical shell for all radii. cylRad[ns] is the cylindrical radius. output=[arr1[ns], …, arrN[ns]] contains the z-integrated output arrays.

Return type:

tuple

class magic.Butterfly(file=None, step=1, iplot=True, rad=0.8, lastvar=None, nvar='all', levels=20, cm='RdYlBu_r', precision=<class 'numpy.float32'>, cut=0.8)[source]

This class can be used to display the time evolution of the magnetic field for various latitudes (i.e. the well-known butterfly diagrams). These diagrams are usually constructed using MagIC’s movie files: either radial cuts (like Br_CMB_mov.TAG) or azimuthal-average (like AB_mov.TAG).

>>> # Read Br_CMB_mov.ccondAnelN3MagRa2e7Pm2ggg
>>> t1 = Butterfly(file='Br_CMB_mov.ccondAnelN3MagRa2e7Pm2ggg', step=1,
                   iplot=False)
>>> # Plot it
>>> t1.plot(levels=33, cm='seismic', cut=0.6)
__add__(new)[source]

Overload of the addition operator

>>> # Read 2 files
>>> b1 = Butterfly(file='AB_mov.test1', iplot=False)
>>> b2 = Butterfly(file='AB_mov.test2', iplot=False)
>>> # Stack them and display the whole thing
>>> b = b1+b2
>>> b.plot(levels=33, contour=True, cut=0.8, cm='seismic')
__init__(file=None, step=1, iplot=True, rad=0.8, lastvar=None, nvar='all', levels=20, cm='RdYlBu_r', precision=<class 'numpy.float32'>, cut=0.8)[source]
Parameters:
  • file (str) – when specified, the constructor reads this file, otherwise a list with the possible options is displayed

  • rad (float) – radial level (normalised to the outer boundary radius)

  • iplot (bool) – display/hide the plots (default is True)

  • nvar (int) – the number of time steps (lines) of the movie file one wants to plot starting from the last line

  • lastvar (int) – the number of the last time step to be read

  • step (int) – the stepping between two lines

  • levels (int) – the number of contour levels

  • cm (str) – the name of the color map

  • cut (float) – adjust the contour extrema to max(abs(data))*cut

  • precision (bool) – precision of the input file, np.float32 for single precision, np.float64 for double precision

__weakref__

list of weak references to the object (if defined)

fourier2D(renorm=False)[source]

This function allows to conduct some basic Fourier analysis on the data. It displays two figures: the first one is a contour levels in the (Frequency, Latitude) plane, the second one is integrated over latitudes (thus a simple, power vs Frequency plot)

>>> # Load the data without plotting
>>> b1 = Butterfly(file='AB_mov.test1', iplot=False)
>>> # Fourier analysis
>>> b1.fourier2D()
Parameters:

renorm (bool) – when set to True, it rebins the time series in case of irregularly spaced data

plot(levels=12, contour=False, renorm=False, cut=0.5, mesh=3, cm='RdYlBu_R')[source]

Plotting function

Parameters:
  • cm (str) – name of the colormap

  • levels (int) – the number of contour levels (only used when iplot=True and contour=True)

  • contour (bool) – when set to True, display contour levels (pylab.contourf), when set to False, display an image (pylab.imshow)

  • renorm (bool) – when set to True, it re-bins the time series in case of irregularly time-spaced data

  • mesh (int) – when renorm=True, factor of regriding: NewTime = mesh*OldTime

  • cut (float) – adjust the contour extrema to max(abs(data))*cut