feathern
commented
1 week ago Thanks for bringing this to our attention. Since it's the middle of the summer, a lot of us have been traveling between conferences and vacation -- which is to say sorry for the delayed response. I'm back now and will have a look at this in the coming week. I'm also going to tag @illorenzo7 on this since he's been making some modifications to the custom-reference-state interface during the past year, and he may see or think of something obvious that doesn't jump out at me right away. -- Nick
Ayesha714
I am running a case with custom_reference_type = 4 using this Python file. However, my code is not reading this file correctly (both constants and functions are being read incorrectly). Could you please guide me on what might be the issue?
This code provides an example for using a custom non-dimensionalization of Rayleigh in Boussinesq dynamo mode.
The non-dimensionalization used here is based on the Rayleigh default viscous diffusion scaling for convection, thus providing a check on the custom reference state.
The parameters correspond to a case listed in Table 2 of Soderlund et al.: "The influence of magnetic fields in planetary dynamo models", Earth Planet. Sci. Lett., v.333-334, p.9-20 (2012)
The numbers referenced below for the various functions and constants refer to equation (5) in "Rayleigh_Output_Variables.pdf". Please refer to that document for further details.
Requirements: (1) "rayleigh_diagnostics.py" ; and (2) "reference_tools.py"
'''
import numpy as np
import os, sys sys.path.insert(0, os.path.abspath('../../'))
import post_processing.reference_tools as rt
name of output file containing custom reference data
filename = 'custom_ref_viscous.dat'our issue content here.
Non-dimensional input parameters
Ra = 1.12e5 # Rayleigh number Pr = 1.0 # Prandtl number Ek = 2.0e-3 # Ekman number Pm = 5.0 # Magnetic Prandtl number beta = 0.4 # Aspect ratio = r_inner/r_outer gravity_power = 1.0 # power law variation of gravitational acceleration
Create radial grid
Numer of radial grid points for radial functions (f_1, f_2, etc.)
Make large enough for accurate interpolation onto Chebyshev grid
nr = 2000
non-dimensional r_inner
ri = beta/(1-beta)
non-dimensional r_outer
ro=1.0/(1-beta)
non-dimensional radial grid
radius=np.linspace(ri,ro,nr)
Define the reference state functions and constants
ones = np.ones(nr,dtype='float64') zeros = np.zeros(nr,dtype='float64')
the function list below is default for Boussinesq
f_1 = ones f_2 = (radius/radius[nr-1])**gravity_power f_3 = ones f_4 = ones f_5 = ones f_6 = zeros f_7 = ones f_8 = zeros f_9 = zeros f_10 = zeros f_11 = zeros f_12 = zeros f_13 = zeros
c_1 = 2.0/Ek # Coriolis force c_2 = Ra/Pr # Buoyancy force c_3 = 1.0/Ek # Pressure gradient c_4 = 1.0/(Ek*Pm) # Lorentz force c_5 = 1.0 # Viscous force c_6 = 1.0/Pr # Thermal diffusion c_7 = 1.0/Pm # Ohmic diffusion c_8 = 0.0 # Viscous heating c_9 = 0.0 # Ohmic heating c_10 = 0.0 # Internal heating
Set all of the functions and constants
my_ref = rt.equation_coefficients(radius)
Set functions here
my_ref.set_function(f_1, 1) my_ref.set_function(f_2, 2) my_ref.set_function(f_3, 3) my_ref.set_function(f_4, 4) my_ref.set_function(f_5, 5) my_ref.set_function(f_6, 6) my_ref.set_function(f_7, 7) my_ref.set_function(f_8, 8) my_ref.set_function(f_9, 9) my_ref.set_function(f_10, 10) my_ref.set_function(f_11, 11) my_ref.set_function(f_12, 12) my_ref.set_function(f_13, 13)
Set constants here
my_ref.set_constant(c_1, 1) my_ref.set_constant(c_2, 2) my_ref.set_constant(c_3, 3) my_ref.set_constant(c_4, 4) my_ref.set_constant(c_5, 5) my_ref.set_constant(c_6, 6) my_ref.set_constant(c_7, 7) my_ref.set_constant(c_8, 8) my_ref.set_constant(c_9, 9) my_ref.set_constant(c_10, 10)
my_ref.write(filename)
print('Custom reference file', filename, 'was written successfully.')