gafusion / omas

Ordered Multidimensional Array Structure
http://gafusion.github.io/omas
MIT License
30 stars 15 forks source link

Possible interpolation issue #148

Closed orso82 closed 3 years ago

orso82 commented 3 years ago

mention @jmcclena

Noticed some difference in density and density_thermal (no density_fast) after calling physics_core_profiles_consistent in an omas_environment where the rho coordinate was coming from a 129 grid equilibrium.

This should reproduce the issue:

from omfit_classes.utils_fusion import Hmode_profiles
rho_eq = np.linspace(0,1,129)**.5
psi_eq = np.linspace(0,1,129)
rho = np.linspace(0,1,201)

ods = ODS()
coordsio = {'core_profiles.profiles_1d.%d.grid.rho_tor_norm' % 0: rho_eq}

ods['core_profiles.profiles_1d[0].grid.rho_tor_norm'] = rho
nval = Hmode_profiles(rgrid=129)
with omas_environment(ods, coordsio=coordsio):
    ods['core_profiles.profiles_1d[0].electrons.density'] = nval
    ods.physics_core_profiles_consistent()

image

orso82 commented 3 years ago

Problem is due to use of linear interpolation. Will need to add an option to specify order of the interpolation

from omfit_classes.utils_fusion import Hmode_profiles
rho_eq = np.linspace(0,1,129)**.5
psi_eq = np.linspace(0,1,129)
rho = np.linspace(0,1,201)

ods = ODS()
coordsio = {'core_profiles.profiles_1d.%d.grid.rho_tor_norm' % 0: rho_eq}

ods['core_profiles.profiles_1d[0].grid.rho_tor_norm'] = rho
nval = Hmode_profiles(rgrid=129)
with omas_environment(ods, coordsio=coordsio):
    ods['core_profiles.profiles_1d[0].electrons.density'] = nval #interpolates 129 to 201
    tmp=ods['core_profiles.profiles_1d[0].electrons.density'] #interpolates 201 to 129
    ods['core_profiles.profiles_1d[0].electrons.density_thermal']=tmp #interpolates 129 to 201

plot(ods['core_profiles.profiles_1d[0].electrons.density'])
plot(ods['core_profiles.profiles_1d[0].electrons.density_thermal'])