Exact solution used in MES runs

We would like to MES the operation

$$ \iiint f \text{d}V = \int_0^{Lx} \int_0^{Ly} \int_0^{2\pi} f \rho \text{d}\theta \text{d}y \text{d}\rho $$

Using cylindrical geometry.



In [1]:

    
%matplotlib notebook
import numpy as np

from sympy import init_printing
from sympy import S
from sympy import sin, cos, tanh, exp, pi, sqrt
from sympy import integrate

from boutdata.mms import x, y, z, t

import os, sys
# If we add to sys.path, then it must be an absolute path
common_dir = os.path.abspath('./../../../../common')
# Sys path is a list of system paths
sys.path.append(common_dir)
from CELMAPy.MES import get_metric, BOUT_print

init_printing()

Initialize



In [2]:

    
folder = '../sumSines/'
metric = get_metric()

Define the variables



In [3]:

    
# Initialization
the_vars = {}

Define the function to volume integrate

NOTE:

These do not need to be fulfilled in order to get convergence

z must be periodic
The field $f(\rho, \theta)$ must be of class infinity in $z=0$ and $z=2\pi$
The field $f(\rho, \theta)$ must be continuous in the $\rho$ direction with $f(\rho, \theta + \pi)$

But this needs to be fulfilled:

The field $f(\rho, \theta)$ must be single valued when $\rho\to0$
Eventual BC in $\rho$ must be satisfied



In [4]:

    
# We need Lx and Ly
from boututils.options import BOUTOptions
myOpts = BOUTOptions(folder)
Lx = eval(myOpts.geom['Lx'])
Ly = eval(myOpts.geom['Ly'])



In [5]:

    
# z sin function

# NOTE: The function is not continuous over origo

s = 2
c = pi
w = pi/2
the_vars['f'] = sin(z + 2.0) + sin(2.0*z + 2.1) + sin(3.0*z)

Calculating the solution



In [6]:

    
rhoInt        =  integrate(the_vars['f']*x, (x, 0 ,Lx))
rhoYInt       =  integrate(rhoInt         , (y, 0 ,Ly))
the_vars['S'] = (integrate(rhoYInt        , (z, 0, 2*np.pi))).evalf()

Print the variables in BOUT++ format



In [7]:

    
BOUT_print(the_vars, rational=False)