IO: RMS force balance, torsional oscillations, misc

RMS.f90

Description

This module contains the calculation of the RMS force balance and induction terms.

Quick access

Variables:

adv2hint, advp2, advp2lm, advrmsl, advrmslnr, advt2, advt2lm, arc2hint, arcmag2hint, arcmagrmsl, arcmagrmslnr, arcrmsl, arcrmslnr, buo_temp2hint, buo_temprmsl, buo_temprmslnr, buo_xi2hint, buo_xirmsl, buo_xirmslnr, cfp2, cfp2lm, cft2, cft2lm, cia2hint, ciarmsl, ciarmslnr, clf2hint, clfrmsl, clfrmslnr, cor2hint, corrmsl, corrmslnr, difpol2hint, difpollmr, difrmsl, difrmslnr, diftor2hint, dpdpc, dpdtc, dpkindrc, dpkindrlm, dr_facc, dtbpol2hint, dtbpollmr, dtbrms_file, dtbtor2hint, dtvp, dtvplm, dtvr, dtvrlm, dtvrms_file, dtvt, dtvtlm, geo2hint, geormsl, geormslnr, iner2hint, inerrmsl, inerrmslnr, lf2hint, lfp2, lfp2lm, lfrlm, lfrmsl, lfrmslnr, lft2, lft2lm, mag2hint, magrmsl, magrmslnr, n_cheb_maxc, n_dtbrms_file, n_dtvrms_file, n_r_maxc, ncut, pfp2lm, pft2lm, plf2hint, plfrmsl, plfrmslnr, pre2hint, prermsl, prermslnr, rc, vp_old, vr_old, vt_old

Routines:

compute_lm_forces(), dtbrms(), dtvrms(), finalize_rms(), get_force(), get_nl_rms(), init_rnb(), initialize_rms(), transform_to_grid_rms(), transform_to_lm_rms(), zerorms()

Needed modules

Variables

  • rms/adv2hint (*,*) [real,private/allocatable]
  • rms/advp2 (*,*) [real,private/allocatable]
  • rms/advp2lm (*) [complex,private/allocatable]
  • rms/advrmsl [mean_sd_type,private]
  • rms/advrmslnr [mean_sd_2d_type,private]
  • rms/advt2 (*,*) [real,private/allocatable]
  • rms/advt2lm (*) [complex,private/allocatable]
  • rms/arc2hint (*,*) [real,private/allocatable]
  • rms/arcmag2hint (*,*) [real,private/allocatable]
  • rms/arcmagrmsl [mean_sd_type,private]
  • rms/arcmagrmslnr [mean_sd_2d_type,private]
  • rms/arcrmsl [mean_sd_type,private]
  • rms/arcrmslnr [mean_sd_2d_type,private]
  • rms/buo_temp2hint (*,*) [real,private/allocatable]
  • rms/buo_temprmsl [mean_sd_type,private]
  • rms/buo_temprmslnr [mean_sd_2d_type,private]
  • rms/buo_xi2hint (*,*) [real,private/allocatable]
  • rms/buo_xirmsl [mean_sd_type,private]
  • rms/buo_xirmslnr [mean_sd_2d_type,private]
  • rms/cfp2 (*,*) [real,private/allocatable]
  • rms/cfp2lm (*) [complex,private/allocatable]
  • rms/cft2 (*,*) [real,private/allocatable]
  • rms/cft2lm (*) [complex,private/allocatable]
  • rms/cia2hint (*,*) [real,private/allocatable]
  • rms/ciarmsl [mean_sd_type,private]
  • rms/ciarmslnr [mean_sd_2d_type,private]
  • rms/clf2hint (*,*) [real,private/allocatable]
  • rms/clfrmsl [mean_sd_type,private]
  • rms/clfrmslnr [mean_sd_2d_type,private]
  • rms/cor2hint (*,*) [real,private/allocatable]
  • rms/corrmsl [mean_sd_type,private]
  • rms/corrmslnr [mean_sd_2d_type,private]
  • rms/difpol2hint (*,*) [real,allocatable/public]
  • rms/difpollmr (*,*) [complex,allocatable/public]
  • rms/difrmsl [mean_sd_type,private]
  • rms/difrmslnr [mean_sd_2d_type,private]
  • rms/diftor2hint (*,*) [real,allocatable/public]
  • rms/dpdpc (*,*) [real,private/allocatable]
  • rms/dpdtc (*,*) [real,private/allocatable]
  • rms/dpkindrc (*,*) [real,private/allocatable]
  • rms/dpkindrlm (*) [complex,private/allocatable]
  • rms/dr_facc (*) [real,allocatable/public]
  • rms/dtbpol2hint (*,*) [real,allocatable/public]
  • rms/dtbpollmr (*,*) [complex,allocatable/public]
  • rms/dtbrms_file [character(len=72),private]
  • rms/dtbtor2hint (*,*) [real,allocatable/public]
  • rms/dtvp (*,*) [real,private/allocatable]
  • rms/dtvplm (*) [complex,private/allocatable]
  • rms/dtvr (*,*) [real,private/allocatable]
  • rms/dtvrlm (*) [complex,private/allocatable]
  • rms/dtvrms_file [character(len=72),private]
  • rms/dtvt (*,*) [real,private/allocatable]
  • rms/dtvtlm (*) [complex,private/allocatable]
  • rms/geo2hint (*,*) [real,private/allocatable]
  • rms/geormsl [mean_sd_type,private]
  • rms/geormslnr [mean_sd_2d_type,private]
  • rms/iner2hint (*,*) [real,private/allocatable]
  • rms/inerrmsl [mean_sd_type,private]
  • rms/inerrmslnr [mean_sd_2d_type,private]
  • rms/lf2hint (*,*) [real,private/allocatable]
  • rms/lfp2 (*,*) [real,private/allocatable]
  • rms/lfp2lm (*) [complex,private/allocatable]
  • rms/lfrlm (*) [complex,private/allocatable]
  • rms/lfrmsl [mean_sd_type,private]
  • rms/lfrmslnr [mean_sd_2d_type,private]
  • rms/lft2 (*,*) [real,private/allocatable]
  • rms/lft2lm (*) [complex,private/allocatable]
  • rms/mag2hint (*,*) [real,private/allocatable]
  • rms/magrmsl [mean_sd_type,private]
  • rms/magrmslnr [mean_sd_2d_type,private]
  • rms/n_cheb_maxc [integer,public]

    Number of Chebyshevs

  • rms/n_dtbrms_file [integer,private]
  • rms/n_dtvrms_file [integer,private]
  • rms/n_r_maxc [integer,public]

    Number of radial points

  • rms/ncut [integer,public]

    Number of points for the cut-off

  • rms/pfp2lm (*) [complex,private/allocatable]
  • rms/pft2lm (*) [complex,private/allocatable]
  • rms/plf2hint (*,*) [real,private/allocatable]
  • rms/plfrmsl [mean_sd_type,private]
  • rms/plfrmslnr [mean_sd_2d_type,private]
  • rms/pre2hint (*,*) [real,private/allocatable]
  • rms/prermsl [mean_sd_type,private]
  • rms/prermslnr [mean_sd_2d_type,private]
  • rms/rc (*) [real,private/allocatable]

    Cut-off radii

  • rms/vp_old (*,*,*) [real,private/allocatable]
  • rms/vr_old (*,*,*) [real,private/allocatable]
  • rms/vt_old (*,*,*) [real,private/allocatable]

Subroutines and functions

subroutine  rms/initialize_rms()

Memory allocation of arrays used in the computation of r.m.s. force balance

Called from:

magic

Call to:

init_rnb(), zerorms()

subroutine  rms/finalize_rms()

Deallocate arrays used for r.m.s. force balance computation

Called from:

magic

subroutine  rms/zerorms()

Zeros integrals that are set in get_td, update_z, update_wp, update_b, dtVrms and dtBrms

Called from:

initialize_rms(), output()

subroutine  rms/init_rnb(r, rcut, rdea, r2, n_r_max2, n_cheb_max2, ns, rscheme_rms)

Prepares the usage of a cut back radial grid where nS points on both boundaries are discarded. The aim actually is to discard boundary effects, but just not considering the boundary grid points does not work when you also want radial derivatives and integrals. For these we use the Chebychev transform which needs are particular number of grid points so that the fast cosine transform can be applied. Therefor more than just 2 points have to be thrown away, which may make sense anyway.

Parameters:
  • r (n_r_max) [real ,in]

  • rcut [real ,in]

  • rdea [real ,in]

  • r2 (*) [real ,out,allocatable]

  • n_r_max2 [integer ,out] :: ‘)

  • n_cheb_max2 [integer ,out]

  • ns [integer ,out] :: ‘)

  • rscheme_rms [real ]

Called from:

initialize_rms()

Call to:

abortrun()

subroutine  rms/get_nl_rms(nr, vr, vt, vp, dvrdr, dvrdt, dvrdp, dvtdr, dvtdp, dvpdr, dvpdp, cvr, advt, advp, lft, lfp, tscheme, lrmscalc)

This subroutine computes the r.m.s. force balance terms which need to be computed on the grid

Parameters:
  • nr [integer ,in] :: radial level

  • vr (*,*) [real ,in]

  • vt (*,*) [real ,in]

  • vp (*,*) [real ,in]

  • dvrdr (*,*) [real ,in]

  • dvrdt (*,*) [real ,in]

  • dvrdp (*,*) [real ,in]

  • dvtdr (*,*) [real ,in]

  • dvtdp (*,*) [real ,in]

  • dvpdr (*,*) [real ,in]

  • dvpdp (*,*) [real ,in]

  • cvr (*,*) [real ,in]

  • advt (*,*) [real ,in]

  • advp (*,*) [real ,in]

  • lft (*,*) [real ,in]

  • lfp (*,*) [real ,in]

  • tscheme [real ] :: time scheme

  • lrmscalc [logical ,in]

Called from:

radialloop()

Call to:

get_openmp_blocks()

subroutine  rms/transform_to_grid_rms(nr, p_rloc)

This subroutine is used to transform the arrays used in r.m.s. force calculations from the spectral space to the physical grid.

Parameters:
  • nr [integer ,in] :: radial level

  • p_rloc (lm_max,1 - nrstart + nrstop) [complex ,inout] :: pressure in LM space

Called from:

transform_to_grid_space()

Call to:

scal_to_grad_spat()

subroutine  rms/transform_to_lm_rms(nr, lfr)

This subroutine is used to transform the arrays used in r.m.s. force calculations from the physical grid to spectral space.

Parameters:
  • nr [integer ,in] :: radial level

  • lfr (*,*) [real ,inout] :: radial component of the Lorentz force

Called from:

transform_to_lm_space()

Call to:

scal_to_sh(), spat_to_sphertor(), spat_to_qst()

subroutine  rms/compute_lm_forces(nr, w_rloc, dw_rloc, ddw_rloc, z_rloc, s_rloc, xi_rloc, p_rloc, dp_rloc, advrlm)

This subroutine finalizes the computation of the r.m.s. spectra once the quantities are back in spectral space.

Parameters:
  • nr [integer ,in]

  • w_rloc (*) [complex ,in]

  • dw_rloc (*) [complex ,in] :: phi-deriv of dw/dr

  • ddw_rloc (*) [complex ,in]

  • z_rloc (*) [complex ,in]

  • s_rloc (*) [complex ,in]

  • xi_rloc (*) [complex ,in]

  • p_rloc (*) [complex ,in]

  • dp_rloc (*) [complex ,in]

  • advrlm (*) [complex ,in]

Called from:

radialloop()

Call to:

hintrms()

subroutine  rms/get_force(force2hint, forcerms, forcermsl, forcermslnr, volc, nrms_sets, timepassed, timenorm, l_stop_time[, forcepol2hint[, forcetor2hint]])

This subroutine is used to compute the contributions of the forces in the Navier-Stokes equation

Parameters:
  • force2hint (1 + l_max,n_r_max) [real ,in]

  • forcerms [real ,out]

  • forcermsl [mean_sd_type ,inout]

  • forcermslnr [mean_sd_2d_type ,inout]

  • volc [real ,in]

  • nrms_sets [integer ,in]

  • timepassed [real ,in]

  • timenorm [real ,in]

  • l_stop_time [logical ,in]

Options:
Called from:

dtvrms()

Call to:

rint_r()

subroutine  rms/dtvrms(time, nrms_sets, timepassed, timenorm, l_stop_time)

This routine calculates and stores the different contributions of the forces entering the Navier-Stokes equation.

Parameters:
  • time [real ,in]

  • nrms_sets [integer ,inout]

  • timepassed [real ,in]

  • timenorm [real ,in]

  • l_stop_time [logical ,in]

Called from:

output()

Call to:

get_dr_rloc(), hint2dpol(), get_force()

subroutine  rms/dtbrms(time)
Parameters:

time [real ,in]

Called from:

output()

Call to:

get_poltorrms(), get_dr_rloc(), hint2dpollm(), rint_r()

RMS_helpers.f90

Description

This module contains several useful subroutines required to compute RMS diagnostics

Quick access

Routines:

get_poltorrms(), hint2dpol(), hint2dpollm(), hintrms()

Needed modules

Variables

Subroutines and functions

subroutine  rms_helpers/get_poltorrms(pol, drpol, tor, llm, ulm, polrms, torrms, polasrms, torasrms, map)

calculates integral PolRms=sqrt( Integral (pol^2 dV) ) calculates integral TorRms=sqrt( Integral (tor^2 dV) ) plus axisymmetric parts. integration in theta,phi by summation of spherical harmonics integration in r by using Chebycheff integrals The mapping map gives the mapping lm to l,m for the input arrays Pol,drPol and Tor Output: PolRms,TorRms,PolAsRms,TorAsRms

Parameters:
  • pol (1 - llm + ulm,n_r_max) [complex ,in] :: Poloidal field Potential

  • drpol (1 - llm + ulm,n_r_max) [complex ,in] :: Radial derivative of Pol

  • tor (1 - llm + ulm,n_r_max) [complex ,in] :: Toroidal field Potential

  • llm [integer ,in,required]

  • ulm [integer ,in,required]

  • polrms [real ,out]

  • torrms [real ,out]

  • polasrms [real ,out]

  • torasrms [real ,out]

  • map [mappings ,in]

Called from:

dtbrms()

Call to:

cc2real(), rint_r()

subroutine  rms_helpers/hint2dpol(dpol, lmstart, lmstop, pol2hint, map)
Parameters:
  • dpol (1 - lmstart + lmstop) [complex ,in] :: Toroidal field Potential

  • lmstart [integer ,in,required]

  • lmstop [integer ,in,required]

  • pol2hint (1 + l_max) [real ,inout]

  • map [mappings ,in]

Called from:

dtvrms()

Call to:

cc2real()

subroutine  rms_helpers/hint2dpollm(dpol, lmstart, lmstop, pol2hint, map)
Parameters:
  • dpol (1 - lmstart + lmstop) [complex ,in]

  • lmstart [integer ,in,required]

  • lmstop [integer ,in,required]

  • pol2hint (1 - lmstart + lmstop) [real ,inout]

  • map [mappings ,in]

Called from:

dtbrms()

Call to:

cc2real()

subroutine  rms_helpers/hintrms(f, nr, lmstart, lmstop, f2hint, map, sphertor)
Parameters:
  • f (*) [complex ,in]

  • nr [integer ,in]

  • lmstart [integer ,in]

  • lmstop [integer ,in]

  • f2hint (1 + l_max) [real ,inout]

  • map [mappings ,in]

  • sphertor [logical ,in]

Called from:

compute_lm_forces()

Call to:

cc2real()

dtB.f90

Description

This module contains magnetic field stretching and advection terms plus a separate omega-effect. It is used for movie output.

Quick access

Variables:

dtb_lmloc_container, dtb_rloc_container, padvlm, padvlm_lmloc, padvlm_rloc, padvlmic, padvlmic_lmloc, pdiflm, pdiflm_lmloc, pdiflmic, pdiflmic_lmloc, pstrlm, pstrlm_lmloc, pstrlm_rloc, tadvlm, tadvlm_lmloc, tadvlm_rloc, tadvlmic, tadvlmic_lmloc, tadvrlm_lmloc, tadvrlm_rloc, tdiflm, tdiflm_lmloc, tdiflmic, tdiflmic_lmloc, tomelm, tomelm_lmloc, tomelm_rloc, tomerlm_lmloc, tomerlm_rloc, tstrlm, tstrlm_lmloc, tstrlm_rloc, tstrrlm_lmloc, tstrrlm_rloc

Routines:

dtb_gather_lo_on_rank0(), finalize_dtb_mod(), get_dh_dtblm(), get_dtblm(), get_dtblmfinish(), initialize_dtb_mod()

Needed modules

Variables

  • dtb_mod/dtb_lmloc_container (*,*,*) [complex,private/target/allocatable]
  • dtb_mod/dtb_rloc_container (*,*,*) [complex,private/target/allocatable]
  • dtb_mod/padvlm (*,*) [complex,allocatable/public]
  • dtb_mod/padvlm_lmloc (*,*) [complex,pointer/public]
  • dtb_mod/padvlm_rloc (*,*) [complex,private/pointer]
  • dtb_mod/padvlmic (*,*) [complex,allocatable/public]
  • dtb_mod/padvlmic_lmloc (*,*) [complex,allocatable/public]
  • dtb_mod/pdiflm (*,*) [complex,allocatable/public]
  • dtb_mod/pdiflm_lmloc (*,*) [complex,allocatable/public]
  • dtb_mod/pdiflmic (*,*) [complex,allocatable/public]
  • dtb_mod/pdiflmic_lmloc (*,*) [complex,allocatable/public]
  • dtb_mod/pstrlm (*,*) [complex,allocatable/public]
  • dtb_mod/pstrlm_lmloc (*,*) [complex,pointer/public]
  • dtb_mod/pstrlm_rloc (*,*) [complex,private/pointer]
  • dtb_mod/tadvlm (*,*) [complex,allocatable/public]
  • dtb_mod/tadvlm_lmloc (*,*) [complex,pointer/public]
  • dtb_mod/tadvlm_rloc (*,*) [complex,private/pointer]
  • dtb_mod/tadvlmic (*,*) [complex,allocatable/public]
  • dtb_mod/tadvlmic_lmloc (*,*) [complex,allocatable/public]
  • dtb_mod/tadvrlm_lmloc (*,*) [complex,pointer/public]
  • dtb_mod/tadvrlm_rloc (*,*) [complex,private/pointer]
  • dtb_mod/tdiflm (*,*) [complex,allocatable/public]
  • dtb_mod/tdiflm_lmloc (*,*) [complex,allocatable/public]
  • dtb_mod/tdiflmic (*,*) [complex,allocatable/public]
  • dtb_mod/tdiflmic_lmloc (*,*) [complex,allocatable/public]
  • dtb_mod/tomelm (*,*) [complex,allocatable/public]
  • dtb_mod/tomelm_lmloc (*,*) [complex,pointer/public]
  • dtb_mod/tomelm_rloc (*,*) [complex,private/pointer]
  • dtb_mod/tomerlm_lmloc (*,*) [complex,pointer/public]
  • dtb_mod/tomerlm_rloc (*,*) [complex,private/pointer]
  • dtb_mod/tstrlm (*,*) [complex,allocatable/public]
  • dtb_mod/tstrlm_lmloc (*,*) [complex,pointer/public]
  • dtb_mod/tstrlm_rloc (*,*) [complex,private/pointer]
  • dtb_mod/tstrrlm_lmloc (*,*) [complex,pointer/public]
  • dtb_mod/tstrrlm_rloc (*,*) [complex,private/pointer]

Subroutines and functions

subroutine  dtb_mod/initialize_dtb_mod()

Memory allocation

Called from:

magic

subroutine  dtb_mod/finalize_dtb_mod()

Memory deallocation

Called from:

magic

subroutine  dtb_mod/dtb_gather_lo_on_rank0()

MPI gather on rank0 for dtBmovie outputs. This routine should really be suppressed once the movie outputs have been improved

Called from:

get_dtblmfinish()

Call to:

gather_all_from_lo_to_rank0()

subroutine  dtb_mod/get_dtblm(nr, vr, vt, vp, br, bt, bp, btvrlm, bpvrlm, brvtlm, brvplm, btvplm, bpvtlm, brvzlm, btvzlm, bpvtbtvpcotlm, bpvtbtvpsn2lm, btvzsn2lm)

This subroutine calculates non-linear products in grid-space for radial level nR.

Parameters:
  • nr [integer ,in]

  • vr (*,*) [real ,in]

  • vt (*,*) [real ,in]

  • vp (*,*) [real ,in]

  • br (*,*) [real ,in]

  • bt (*,*) [real ,in]

  • bp (*,*) [real ,in]

  • btvrlm (*) [complex ,out]

  • bpvrlm (*) [complex ,out]

  • brvtlm (*) [complex ,out]

  • brvplm (*) [complex ,out]

  • btvplm (*) [complex ,out]

  • bpvtlm (*) [complex ,out]

  • brvzlm (*) [complex ,out]

  • btvzlm (*) [complex ,out]

  • bpvtbtvpcotlm (*) [complex ,out]

  • bpvtbtvpsn2lm (*) [complex ,out]

  • btvzsn2lm (*) [complex ,out]

Called from:

radialloop()

Call to:

spat_to_sphertor(), scal_to_sh()

subroutine  dtb_mod/get_dh_dtblm(nr, btvrlm, bpvrlm, brvtlm, brvplm, btvplm, bpvtlm, brvzlm, btvzlm, bpvtbtvpcotlm, bpvtbtvpsn2lm)

Purpose of this routine is to calculate theta and phi derivative related terms of the magnetic production and advection terms and store them.

Parameters:
  • nr [integer ,in]

  • btvrlm (*) [complex ,in]

  • bpvrlm (*) [complex ,in]

  • brvtlm (*) [complex ,in]

  • brvplm (*) [complex ,in]

  • btvplm (*) [complex ,in]

  • bpvtlm (*) [complex ,in]

  • brvzlm (*) [complex ,in]

  • btvzlm (*) [complex ,in]

  • bpvtbtvpcotlm (*) [complex ,in]

  • bpvtbtvpsn2lm (*) [complex ,in]

Called from:

radialloop()

subroutine  dtb_mod/get_dtblmfinish(time, n_time_step, omega_ic, b, ddb, aj, dj, ddj, b_ic, db_ic, ddb_ic, aj_ic, dj_ic, ddj_ic, l_frame)

– Input of variables:

Parameters:
Called from:

output()

Call to:

get_openmp_blocks(), rbrspec(), rbpspec(), dtb_gather_lo_on_rank0()

dtB_arrays.f90

Quick access

Types:

dtb_arrays_t

Needed modules

Types

  • type  dtb_arrays_mod/dtb_arrays_t
    Type fields:
    • % bpvrlm (*) [complex ,allocatable]

    • % bpvtbtvpcotlm (*) [complex ,allocatable]

    • % bpvtbtvpsn2lm (*) [complex ,allocatable]

    • % bpvtlm (*) [complex ,allocatable]

    • % brvplm (*) [complex ,allocatable]

    • % brvtlm (*) [complex ,allocatable]

    • % brvzlm (*) [complex ,allocatable]

    • % btvplm (*) [complex ,allocatable]

    • % btvrlm (*) [complex ,allocatable]

    • % btvzlm (*) [complex ,allocatable]

    • % btvzsn2lm (*) [complex ,allocatable]

Variables

Subroutines and functions

subroutine  dtb_arrays_mod/initialize(this)
Parameters:

this [real ]

subroutine  dtb_arrays_mod/finalize(this)
Parameters:

this [real ]

subroutine  dtb_arrays_mod/set_zero(this)
Parameters:

this [real ]

out_dtB_frame.f90

Quick access

Routines:

get_bpol(), get_btor(), get_dtb(), lm2pt(), write_dtb_frame()

Needed modules

Variables

Subroutines and functions

subroutine  out_dtb_frame/write_dtb_frame(n_movie, b, db, aj, dj, b_ic, db_ic, aj_ic, dj_ic)

Controls output of specific movie frames related to magnetic field production and diffusion.

Parameters:
Called from:

write_movie_frame()

Call to:

get_dtb(), get_bpol(), get_btor(), lm2pt(), get_drns_even()

subroutine  out_dtb_frame/get_dtb(dtb, dtblm, dimb1, dimb2, n_r, l_ic)
Parameters:
  • dtb (*) [real ,out] :: Result Field with dim>=n_theta_block

  • dtblm (dimb1,dimb2) [complex ,in]

  • dimb1 [integer ,in,]

  • dimb2 [integer ,in,]

  • n_r [integer ,in] :: No. of radial grid point

  • l_ic [logical ,in] :: =true if inner core field

Called from:

write_dtb_frame()

Call to:

toraxi_to_spat()

subroutine  out_dtb_frame/get_bpol(pollm, dpollm, br, bt, bp, rt, lic)

Purpose of this subroutine is to calculate the components Br, Bt, and Bp of the poloidal magnetic field PolLM (l,m space) at the radial grid point r=rT and the block of theta grid points from n_theta=n_theta_start to n_theta=n_theta_start+n_theta_block-1 and for all phis. For lIC=.true. the inner core field is calculated, to get the IC field for a conducting inner core PolLM has to be the poloidal field at the ICB.

Parameters:
  • pollm (lm_max) [complex ,in] :: field in (l,m)-space for rT

  • dpollm (lm_max) [complex ,in] :: dr field in (l,m)-space for rT

  • br (*,*) [real ,out]

  • bt (*,*) [real ,out]

  • bp (*,*) [real ,out]

  • rt [real ,in] :: radius

  • lic [logical ,in] :: true for inner core, special rDep !

Called from:

write_dtb_frame()

Call to:

torpol_to_spat()

subroutine  out_dtb_frame/get_btor(tlm, bt, bp, rt, lic)

Purpose of this subroutine is to calculate the components Bt and Bp of the toroidal magnetic field Tlm (in l,m space) at the radial grid point r=rT and the block of theta grid points from n_theta=n_theta_start to n_theta=n_theta_start+n_theta_block-1 and for all phis. For lIC=.true. the inner core field is calculated, to get the IC field for a conducting inner core Plm has to be the toroidal field at the ICB.

Parameters:
  • tlm (lm_max) [complex ,in] :: field in (l,m)-space for rT

  • bt (*,*) [real ,out]

  • bp (*,*) [real ,out]

  • rt [real ,in] :: radius

  • lic [logical ,in] :: true for inner core, special rDep !

Called from:

write_dtb_frame()

Call to:

sphtor_to_spat()

subroutine  out_dtb_frame/lm2pt(alm, aij, rt, lic, lrcomp)

Spherical harmonic transform from alm(l,m) to aij(phi,theta) Radial field components are calculated for lrComp=.true. Done within the range [n_theta_min,n_theta_min+n_theta_block-1] Used only for graphic output.

Parameters:
  • alm (*) [complex ,in] :: field in (l,m)-space

  • aij (*,*) [real ,out] :: field in (theta,phi)-space

  • rt [real ,in]

  • lic [logical ,in] :: true for inner core, extra factor !

  • lrcomp [logical ,in] :: true for radial field components

Called from:

write_dtb_frame()

Call to:

scal_to_spat()

TO.f90

Description

This module contains information for TO calculation and output

Quick access

Variables:

bplast, bpsdas_rloc, bpzas_rloc, bpzdas_rloc, bs2as_rloc, bslast, bspas_rloc, bspdas_rloc, bszas_rloc, bzlast, bzpdas_rloc, ddzasl, dzasl, dzastras_rloc, dzcoras_rloc, dzddvpas_rloc, dzddvplmr, dzdvpas_rloc, dzdvplmr, dzlfas_rloc, dzrstras_rloc, dzstras_rloc, v2as_rloc, vas_rloc, zasl

Routines:

finalize_to(), get_pas(), getto(), gettofinish(), gettonext(), initialize_to(), prep_to_axi()

Needed modules

Variables

  • torsional_oscillations/bplast (*,*,*) [real,private/allocatable]
  • torsional_oscillations/bpsdas_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/bpzas_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/bpzdas_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/bs2as_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/bslast (*,*,*) [real,private/allocatable]
  • torsional_oscillations/bspas_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/bspdas_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/bszas_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/bzlast (*,*,*) [real,private/allocatable]
  • torsional_oscillations/bzpdas_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/ddzasl (*,*) [real,allocatable/public]
  • torsional_oscillations/dzasl (*) [real,private/allocatable]
  • torsional_oscillations/dzastras_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/dzcoras_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/dzddvpas_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/dzddvplmr (*,*) [real,private/allocatable]
  • torsional_oscillations/dzdvpas_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/dzdvplmr (*,*) [real,private/allocatable]
  • torsional_oscillations/dzlfas_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/dzrstras_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/dzstras_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/v2as_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/vas_rloc (*,*) [real,allocatable/public]
  • torsional_oscillations/zasl (*) [real,private/allocatable]

Subroutines and functions

subroutine  torsional_oscillations/initialize_to()

Allocate the memory needed

Called from:

magic

subroutine  torsional_oscillations/finalize_to()

Deallocate the memory

Called from:

magic

subroutine  torsional_oscillations/prep_to_axi(z, dz)
Parameters:
  • z (*) [complex ,in]

  • dz (*) [complex ,in]

Called from:

radialloop()

subroutine  torsional_oscillations/getto(vr, vt, vp, cvr, dvpdr, br, bt, bp, cbr, cbt, dzrstrlm, dzastrlm, dzcorlm, dzlflm, dtlast, nr)

This program calculates various axisymmetric linear and nonlinear variables for a radial grid point nR and a theta-block. Input are the fields vr,vt,vp,cvr,dvpdr Output are linear azimuthally averaged field VpAS (flow phi component), VpAS2 (square of flow phi component), V2AS (V*V), and Coriolis force Cor. These are give in (r,theta)-space. Also in (r,theta)-space are azimuthally averaged correlations of non-axisymmetric flow components and the respective squares: Vsp=Vs*Vp,Vzp,Vsz,VspC,VzpC,VszC. These are used to calulcate the respective correlations and Reynolds stress. In addition three output field are given in (lm,r) space: dzRstrLMr,dzAstrLMr,dzCorLM,dzLFLM.

These are used to calculate the total Reynolds stress, advection and viscous stress later. Their calculation retraces the calculations done in the time-stepping part of the code.

Parameters:
  • vr (*,*) [real ,in]

  • vt (*,*) [real ,in]

  • vp (*,*) [real ,in]

  • cvr (*,*) [real ,in]

  • dvpdr (*,*) [real ,in]

  • br (*,*) [real ,in]

  • bt (*,*) [real ,in]

  • bp (*,*) [real ,in]

  • cbr (*,*) [real ,in]

  • cbt (*,*) [real ,in]

  • dzrstrlm (1 + l_max) [real ,out]

  • dzastrlm (1 + l_max) [real ,out]

  • dzcorlm (1 + l_max) [real ,out]

  • dzlflm (1 + l_max) [real ,out]

  • dtlast [real ,in] :: last time step

  • nr [integer ,in] :: radial grid point

Called from:

radialloop()

Call to:

spat_to_sh_axi()

subroutine  torsional_oscillations/gettonext(br, bt, bp, ltonext, ltonext2, dt, dtlast, nr)

Preparing TO calculation by storing flow and magnetic field contribution to build time derivative.

Parameters:
  • br (*,*) [real ,in]

  • bt (*,*) [real ,in]

  • bp (*,*) [real ,in]

  • ltonext [logical ,in]

  • ltonext2 [logical ,in]

  • dt [real ,in]

  • dtlast [real ,in]

  • nr [integer ,in]

Called from:

radialloop()

subroutine  torsional_oscillations/gettofinish(nr, dtlast, dzrstrlm, dzastrlm, dzcorlm, dzlflm)

This program was previously part of getTO(…) It has now been separated to get it out of the theta-block loop.

Parameters:
  • nr [integer ,in]

  • dtlast [real ,in]

  • dzrstrlm (1 + l_max) [real ,in]

  • dzastrlm (1 + l_max) [real ,in]

  • dzcorlm (1 + l_max) [real ,in]

  • dzlflm (1 + l_max) [real ,in]

Called from:

radialloop()

Call to:

get_pas()

subroutine  torsional_oscillations/get_pas(tlm, bp, nr)

Purpose of this subroutine is to calculate the axisymmetric component Bp of an axisymmetric toroidal field Tlm given in spherical harmonic space (1:lmax+1). Unscrambling of theta is also ensured

Parameters:
  • tlm (*) [real ,in] :: field in (l,m)-space for rT

  • bp (*) [real ,out]

  • nr [integer ,in]

Called from:

gettofinish()

Call to:

toraxi_to_spat()

TO_arrays.f90

Quick access

Types:

to_arrays_t

Needed modules

Types

  • type  to_arrays_mod/to_arrays_t
    Type fields:
    • % dzastrlm (*) [real ,allocatable]

    • % dzcorlm (*) [real ,allocatable]

    • % dzlflm (*) [real ,allocatable]

    • % dzrstrlm (*) [real ,allocatable]

Variables

Subroutines and functions

subroutine  to_arrays_mod/initialize(this)
Parameters:

this [real ]

subroutine  to_arrays_mod/finalize(this)
Parameters:

this [real ]

subroutine  to_arrays_mod/set_zero(this)
Parameters:

this [real ]

out_TO.f90

Description

This module handles the writing of TO-related outputs: zonal force balance and z-integrated terms. This is a re-implementation of the spectral method used up to MagIC 5.10, which formerly relies on calculation of Plm on the cylindrical grid. This was quite costly and not portable on large truncations. As such, a local 4th order method is preferred here.

Quick access

Variables:

bpsdas, bpzas, bpzdas, bs2as, bspas, bspdas, bszas, bzpdas, cyl, dzastras, dzcoras, dzddvpas, dzdvpas, dzlfas, dzrstras, dzstras, h, movfile, n_nhs_file, n_s_max, n_s_otc, n_shs_file, n_to_file, n_tomov_file, ntomovsets, oh, tofile, v2as, vas, volcyl_oc

Routines:

cylmean(), finalize_outto_mod(), gather_from_rloc_to_rank0(), get_dds(), get_ds(), initialize_outto_mod(), interp_theta(), outto()

Needed modules

Variables

  • outto_mod/bpsdas (*,*) [real,private/allocatable]
  • outto_mod/bpzas (*,*) [real,private/allocatable]
  • outto_mod/bpzdas (*,*) [real,private/allocatable]
  • outto_mod/bs2as (*,*) [real,private/allocatable]
  • outto_mod/bspas (*,*) [real,private/allocatable]
  • outto_mod/bspdas (*,*) [real,private/allocatable]
  • outto_mod/bszas (*,*) [real,private/allocatable]
  • outto_mod/bzpdas (*,*) [real,private/allocatable]
  • outto_mod/cyl (*) [real,private/allocatable]

    Cylindrical grid

  • outto_mod/dzastras (*,*) [real,private/allocatable]
  • outto_mod/dzcoras (*,*) [real,private/allocatable]
  • outto_mod/dzddvpas (*,*) [real,private/allocatable]
  • outto_mod/dzdvpas (*,*) [real,private/allocatable]
  • outto_mod/dzlfas (*,*) [real,private/allocatable]
  • outto_mod/dzrstras (*,*) [real,private/allocatable]
  • outto_mod/dzstras (*,*) [real,private/allocatable]
  • outto_mod/h (*) [real,private/allocatable]

    height

  • outto_mod/movfile [character(len=64),private]
  • outto_mod/n_nhs_file [integer,private]
  • outto_mod/n_s_max [integer,private]
  • outto_mod/n_s_otc [integer,private]
  • outto_mod/n_shs_file [integer,private]
  • outto_mod/n_to_file [integer,private]
  • outto_mod/n_tomov_file [integer,private]
  • outto_mod/ntomovsets [integer,private]

    Number of TO_mov frames

  • outto_mod/oh (*) [real,private/allocatable]

    1/h

  • outto_mod/tofile [character(len=64),private]
  • outto_mod/v2as (*,*) [real,private/allocatable]
  • outto_mod/vas (*,*) [real,private/allocatable]
  • outto_mod/volcyl_oc [real,private]

Subroutines and functions

subroutine  outto_mod/initialize_outto_mod()

Memory allocation of arrays needed for TO outputs

Called from:

magic

Call to:

simps()

subroutine  outto_mod/finalize_outto_mod()

Memory de-allocation of arrays needed for TO outputs

Called from:

magic

subroutine  outto_mod/outto(time, n_time_step, ekin, ekintas, ltomov)

Output of axisymmetric zonal flow, its relative strength, its time variation, and all forces acting on it.

Parameters:
  • time [real ,in] :: time

  • n_time_step [integer ,in] :: Iteration number

  • ekin [real ,in] :: Kinetic energy

  • ekintas [real ,in] :: Toroidal axisymmetric energy

  • ltomov [logical ,in] :: Do we need to store the movie files as well

Called from:

output()

Call to:

gather_from_rloc_to_rank0(), cylmean(), interp_theta(), get_ds(), get_dds(), simps(), logwrite()

subroutine  outto_mod/cylmean(dat, datn, dats)

This routine computes the z-average inside and outside TC

Parameters:
  • dat (*,*) [real ,in] :: input data

  • datn (n_s_max) [real ,out] :: z-average outside T.C. + N.H.

  • dats (n_s_max) [real ,out] :: z-average outside T.C. + S.H.

Called from:

outgeos(), write_geos_frame(), outto()

Call to:

cylmean_otc(), cylmean_itc()

subroutine  outto_mod/interp_theta(a, ac, rr, cyl, theta)

This routine computes the interpolation of a value at the surface of a spherical shell onto the cylindrical grid. This is only a theta interpolation using the cylindrical theta’s.

Parameters:
  • a (*) [real ,in] :: Field at the outer radius

  • ac (2) [real ,out] :: Surface values for NH and SH as a function of s

  • rr [real ,in] :: Radius at which we compute the extrapolation (can be ri or ro)

  • cyl [real ,in] :: Cylindrical radius

  • theta (*) [real ,in] :: Colatitude

Called from:

outto()

subroutine  outto_mod/get_ds(arr, darr, cyl)

This subroutine is used to compute the 4th order accurate first s-derivative on the regularly spaced grid

Parameters:
  • arr (*) [real ,in] :: Array to be differentiated

  • darr (*) [real ,out] :: s-derivative of the input array

  • cyl (*) [real ,in] :: Cylindrical grid

Called from:

outto()

subroutine  outto_mod/get_dds(arr, ddarr, cyl)

This subroutine is used to compute the 4th order accurate 2nd s-derivative on the regularly spaced grid

https://bellaard.com/tools/Finite%20difference%20coefficient%20calculator/

Parameters:
  • arr (*) [real ,in] :: Array to be differentiated

  • ddarr (*) [real ,out] :: s-derivative of the input array

  • cyl (*) [real ,in] :: Cylindrical grid

Called from:

outto()

subroutine  outto_mod/gather_from_rloc_to_rank0(arr_rloc, arr)

This subroutine gathers the r-distributed array

Parameters:
Called from:

outto()

radial_spectra.f90

Quick access

Variables:

filehandle

Routines:

rbpspec(), rbrspec()

Needed modules

Variables

  • radial_spectra/filehandle [integer,private]

Subroutines and functions

subroutine  radial_spectra/rbrspec(time, pol, polic, fileroot, lic, map)
Parameters:
  • time [real ,in]

  • pol (1 - llm + ulm,n_r_max) [complex ,in]

  • polic (1 - llm + ulm,n_r_ic_max) [complex ,in]

  • fileroot [character(len=*),in]

  • lic [logical ,in]

  • map [mappings ,in]

Called from:

get_dtblmfinish(), output()

Call to:

cc2real()

subroutine  radial_spectra/rbpspec(time, tor, toric, fileroot, lic, map)

Called from rank0, map gives the lm order of Tor and TorIC

Parameters:
  • time [real ,in]

  • tor (1 - llm + ulm,n_r_max) [complex ,in]

  • toric (1 - llm + ulm,n_r_ic_max) [complex ,in]

  • fileroot [character(len=*),in]

  • lic [logical ,in]

  • map [mappings ,in]

Called from:

get_dtblmfinish(), output()

Call to:

cc2real()

outGeos.f90

Description

This module is used to compute z-integrated diagnostics such as the degree of geostrophy or the separation of energies between inside and outside the tangent cylinder. This makes use of a local Simpson’s method. This also handles the computation of the z-average profile of rotation when a Couette flow setup is used

Quick access

Variables:

geos_file, n_geos_file, npstart, npstop, phi_balance, up_ploc, up_rloc, us_ploc, us_rloc, uz_ploc, uz_rloc, vol_otc, wz_ploc, wz_rloc

Routines:

calcgeos(), finalize_geos(), gather_ploc(), initialize_geos(), outgeos(), outomega(), transp_r2phi(), write_geos_frame()

Needed modules

Variables

  • geos/cyl (*) [real,allocatable/public]

    Cylindrical grid

  • geos/geos_file [character(len=72),private]

    file name

  • geos/h (*) [real,private/allocatable]

    h(s)

  • geos/n_geos_file [integer,private]

    file unit for geos.TAG

  • geos/n_s_max [integer,public]

    Number of cylindrical points

  • geos/n_s_otc [integer,private]

    Index for last point outside TC

  • geos/npstart [integer,private]

    Starting nPhi index when MPI distributed

  • geos/npstop [integer,private]

    Stoping nPhi index when MPI distributed

  • geos/phi_balance (*) [load,private/allocatable]

    phi-distributed balance

  • geos/up_ploc (*,*,*) [real,private/allocatable]
  • geos/up_rloc (*,*,*) [real,private/allocatable]
  • geos/us_ploc (*,*,*) [real,private/allocatable]
  • geos/us_rloc (*,*,*) [real,private/allocatable]
  • geos/uz_ploc (*,*,*) [real,private/allocatable]
  • geos/uz_rloc (*,*,*) [real,private/allocatable]
  • geos/vol_otc [real,private]

    volume outside tangent cylinder

  • geos/wz_ploc (*,*,*) [real,private/allocatable]
  • geos/wz_rloc (*,*,*) [real,private/allocatable]

Subroutines and functions

subroutine  geos/initialize_geos(l_geos, l_sric, l_geosmovie)

Memory allocation and definition of the cylindrical grid

Parameters:
  • l_geos [logical ,in] :: Do we need the geos outputs

  • l_sric [logical ,in] :: Is the inner core rotating

  • l_geosmovie [logical ,in] :: Geos movie

Called from:

magic

Call to:

getblocks(), simps()

subroutine  geos/finalize_geos(l_geos, l_sric, l_geosmovie)

Memory deallocation

Parameters:
  • l_geos [logical ,in] :: Do we need the geos outputs

  • l_sric [logical ,in] :: Is the inner core rotating?

  • l_geosmovie [logical ,in] :: Do we have geos movies?

Called from:

magic

subroutine  geos/calcgeos(vr, vt, vp, cvr, dvrdp, dvpdr, nr)

This routine computes the term needed for geos.TAG outputs in physical space.

Parameters:
  • vr (*,*) [real ,in]

  • vt (*,*) [real ,in]

  • vp (*,*) [real ,in]

  • cvr (*,*) [real ,in]

  • dvrdp (*,*) [real ,in]

  • dvpdr (*,*) [real ,in]

  • nr [integer ,in] :: Radial grid point

Called from:

radialloop()

subroutine  geos/outgeos(time, geos, geosa, geosz, geosm, geosnap, ekin)

This routine handles the output of geos.TAG

Parameters:
  • time [real ,in]

  • geos [real ,out]

  • geosa [real ,out]

  • geosz [real ,out]

  • geosm [real ,out]

  • geosnap [real ,out]

  • ekin [real ,out]

Called from:

output()

Call to:

transp_r2phi(), cylmean(), simps(), cylmean_otc()

subroutine  geos/write_geos_frame(n_movie)

This subroutine handles the computation and the writing of geos movie files.

Parameters:

n_movie [integer ,in] :: The index of the movie in list of movies

Called from:

write_movie_frame()

Call to:

transp_r2phi(), cylmean(), gather_ploc()

subroutine  geos/outomega(z, omega_ic)

Output of axisymmetric zonal flow omega(s) into field omega.TAG, where s is the cylindrical radius. This is done for the southern and norther hemispheres at z=+-(r_icb+0.5)

Parameters:
  • z (1 - llm + ulm,n_r_max) [complex ,in] :: Toroidal potential

  • omega_ic [real ,in] :: Rotation rate of the inner core

Called from:

output()

Call to:

toraxi_to_spat(), cylmean_otc(), cylmean_itc()

subroutine  geos/gather_ploc(arr_ploc, arr_full)

This subroutine gathers and transpose a phi-distributed array on rank0.

Parameters:
Called from:

write_geos_frame()

subroutine  geos/transp_r2phi(arr_rloc, arr_ploc)

This subroutine is used to compute a MPI transpose between a R-distributed array and a Phi-distributed array

Parameters:
Called from:

outgeos(), write_geos_frame()

subroutine  geos/cylmean(dat, dat_otc, dat_itc_n, dat_itc_s)

This routine computes the z-average inside and outside TC

Parameters:
  • dat (*,*) [real ,in]

  • dat_otc (n_s_max) [real ,out]

  • dat_itc_n (n_s_max) [real ,out]

  • dat_itc_s (n_s_max) [real ,out]

Call to:

cylmean_otc(), cylmean_itc()

probes.f90

Description

Module for artificial sensors to compare time series of physical data with experiments. Probes are located in a radially symmetrical way on a radial surface given by r_probe (in terms of r_cmb), theta_probe in degrees between 0 and 90 and n_phi_probes denoting the number of probes in phi. Probes will be located at ‘n_phi_probes’ points at two equatorially symmetric latitudes - theta_probe and (180 - theta_probe) on r = r_probe.

version 1.0: Works only for v_phi, for now. Will be extended for other data later.

Quick access

Variables:

n_phi_probes, n_probebr, n_probebt, n_probevp, n_theta_usr, probe_filebr, probe_filebt, probe_filevp, r_probe, rad_rank, rad_usr, theta_probe

Routines:

finalize_probes(), initialize_probes(), probe_out()

Needed modules

Variables

  • probe_mod/n_phi_probes [integer,public]

    number of probes in phi - symmetrically distributed

  • probe_mod/n_probebr [integer,private]
  • probe_mod/n_probebt [integer,private]
  • probe_mod/n_probevp [integer,private]
  • probe_mod/n_theta_usr [integer,private]
  • probe_mod/probe_filebr [character(len=72),private]
  • probe_mod/probe_filebt [character(len=72),private]
  • probe_mod/probe_filevp [character(len=72),private]
  • probe_mod/r_probe [real,public]
  • probe_mod/rad_rank [integer,private]
  • probe_mod/rad_usr [integer,private]
  • probe_mod/theta_probe [real,public]

    probe locations, r_probe in terms of r_cmb and theta in degrees

Subroutines and functions

subroutine  probe_mod/initialize_probes()
Called from:

magic

subroutine  probe_mod/finalize_probes()
Called from:

magic

subroutine  probe_mod/probe_out(time, n_r, vp, br, bt)
Parameters:
  • time [real ,in] :: Time

  • n_r [integer ,in] :: radial grod point no.

  • vp (*,*) [real ,in]

  • br (*,*) [real ,in]

  • bt (*,*) [real ,in]

Called from:

radialloop()