# Path: sucemfem/Consts.py
# Z0 = mu0*c0                              
# 
# Path: sucemfem/Utilities/MeshGenerators.py
# def get_centred_cube(domain_size, max_edge_len, centred_element_coordinate=None):
#     """ Generate a cube mesh centred arount [0,0,0]
# 
#     optionally translates mesh slightly such that
#     centred_element_coordinate is at the centre of an element.
#     """
#     domain_subdivisions = N.array(N.ceil(N.sqrt(2)*domain_size/max_edge_len), N.uint)
#     mesh = dolfin.UnitCube(*domain_subdivisions)
#     # Transform mesh to correct dimensions
#     mesh.coordinates()[:] *= domain_size
#     mesh.coordinates()[:] -= domain_size/2
#     if centred_element_coordinate is not None:
#         ## Translate mesh slightly so that source coordinate lies at
#         ## centroid of an element
#         centred_element_point = dolfin.Point(*centred_element_coordinate)
# #        source_elnos = mesh.all_intersected_entities(centred_element_point)
#         source_elnos = mesh.intersected_cells(centred_element_point)
#         closest_elno = source_elnos[(N.argmin(
#             [centred_element_point.distance(dolfin.Cell(mesh, i).midpoint())
#              for i in source_elnos]))]
#         centre_pt = dolfin.Cell(mesh, closest_elno).midpoint()
#         centre_coord = N.array([centre_pt.x(), centre_pt.y(), centre_pt.z()])
#         # The mesh intersect operator caches certain data
#         # structures. We need to clear them if the mesh coordinates
#         # are changed after calling any mesh intersection methods.
#         mesh.intersection_operator().clear()
#         mesh.coordinates()[:] -= centre_coord
#     return mesh
# 
# Path: sucemfem/PostProcessing/surface_ntff.py
# class SurfaceNTFFForms(object):
# class NTFF(object):
#     def __init__(self, function_space):
#     def set_dofs(self, dofs):
#     def get_N_form(self):
#     def get_L_form(self):
#     def assemble_N(self):
#     def assemble_L(self):
#     def set_parms(self, theta, phi, k0):
#     def __init__(self, function_space):
#     def set_dofs(self, dofs):
#     def set_frequency(self, frequency):
#     def calc_pt(self, theta_deg, phi_deg):
#         N = self.get_N_form()
#         L = self.get_L_form()
# 
# Path: sucemfem/PostProcessing/variational_ntff.py
# class TransformTestingExpression(object):
# class NTFF(object):
#     def __init__(self):
#     def set_parms(self, rhat, ahat, k0):
#     def _set_expr_parms(self, expr, rhat, k0, constfac):
#     def get_expression(self):
#     def __init__(self, function_space, testing_space=None):
#     def set_k0(self, k0):
#     def set_frequency(self, frequency):
#     def set_dofs(self, dofs):
#     def calc_pt(self, theta_deg, phi_deg):
#     def calc_pt_E_H(self, theta_deg, phi_deg):
#     def calc_ff_func(self, rhat, ahat):


import pickle
import numpy as N
import dolfin
import sys
import sucemfem.Utilities.Optimization
import pylab 
from sucemfem.Consts import Z0, c0
from sucemfem.Utilities.MeshGenerators import get_centred_cube
from sucemfem.PostProcessing import surface_ntff, variational_ntff
# Authors:
# Neilen Marais <nmarais@gmail.com>
from __future__ import division


sys.path.append('../../')
#from FenicsCode.Consts import Z0, c0
#from FenicsCode.Utilities.MeshGenerators import get_centred_cube
#import FenicsCode.Utilities.Optimization
#from FenicsCode.PostProcessing import surface_ntff, variational_ntff


# Enable dolfin's form optimizations
sucemfem.Utilities.Optimization.set_dolfin_optimisation()
#FenicsCode.Utilities.Optimization.set_dolfin_optimisation()

# Near-field of an infintesimal dipole
#fname = 'data/dofs-2-0.299792458-0.0166551365556-0.0749481145-1000.pickle'
fname = 'data/dofs-2-0.899377374-0.0499654096667-0.6895226534-1000.pickle'
theta_deg = N.linspace(0, 180, 181)
no_ff_pts = len(theta_deg)
phi_deg = N.zeros(no_ff_pts)
data = pickle.load(open(fname))
lam = c0/data['freq']

# Path: sucemfem/Consts.py
# Z0 = mu0*c0                              
# 
# Path: sucemfem/Utilities/MeshGenerators.py
# def get_centred_cube(domain_size, max_edge_len, centred_element_coordinate=None):
#     """ Generate a cube mesh centred arount [0,0,0]
# 
#     optionally translates mesh slightly such that
#     centred_element_coordinate is at the centre of an element.
#     """
#     domain_subdivisions = N.array(N.ceil(N.sqrt(2)*domain_size/max_edge_len), N.uint)
#     mesh = dolfin.UnitCube(*domain_subdivisions)
#     # Transform mesh to correct dimensions
#     mesh.coordinates()[:] *= domain_size
#     mesh.coordinates()[:] -= domain_size/2
#     if centred_element_coordinate is not None:
#         ## Translate mesh slightly so that source coordinate lies at
#         ## centroid of an element
#         centred_element_point = dolfin.Point(*centred_element_coordinate)
# #        source_elnos = mesh.all_intersected_entities(centred_element_point)
#         source_elnos = mesh.intersected_cells(centred_element_point)
#         closest_elno = source_elnos[(N.argmin(
#             [centred_element_point.distance(dolfin.Cell(mesh, i).midpoint())
#              for i in source_elnos]))]
#         centre_pt = dolfin.Cell(mesh, closest_elno).midpoint()
#         centre_coord = N.array([centre_pt.x(), centre_pt.y(), centre_pt.z()])
#         # The mesh intersect operator caches certain data
#         # structures. We need to clear them if the mesh coordinates
#         # are changed after calling any mesh intersection methods.
#         mesh.intersection_operator().clear()
#         mesh.coordinates()[:] -= centre_coord
#     return mesh
# 
# Path: sucemfem/PostProcessing/surface_ntff.py
# class SurfaceNTFFForms(object):
# class NTFF(object):
#     def __init__(self, function_space):
#     def set_dofs(self, dofs):
#     def get_N_form(self):
#     def get_L_form(self):
#     def assemble_N(self):
#     def assemble_L(self):
#     def set_parms(self, theta, phi, k0):
#     def __init__(self, function_space):
#     def set_dofs(self, dofs):
#     def set_frequency(self, frequency):
#     def calc_pt(self, theta_deg, phi_deg):
#         N = self.get_N_form()
#         L = self.get_L_form()
# 
# Path: sucemfem/PostProcessing/variational_ntff.py
# class TransformTestingExpression(object):
# class NTFF(object):
#     def __init__(self):
#     def set_parms(self, rhat, ahat, k0):
#     def _set_expr_parms(self, expr, rhat, k0, constfac):
#     def get_expression(self):
#     def __init__(self, function_space, testing_space=None):
#     def set_k0(self, k0):
#     def set_frequency(self, frequency):
#     def set_dofs(self, dofs):
#     def calc_pt(self, theta_deg, phi_deg):
#     def calc_pt_E_H(self, theta_deg, phi_deg):
#     def calc_ff_func(self, rhat, ahat):


# Path: examples/filament_source/cube_ntff.py
import pickle
import numpy as N
import dolfin
import sys
import sucemfem.Utilities.Optimization
import pylab 
from sucemfem.Consts import Z0, c0
from sucemfem.Utilities.MeshGenerators import get_centred_cube
from sucemfem.PostProcessing import surface_ntff, variational_ntff

# Authors:
# Neilen Marais <nmarais@gmail.com>
from __future__ import division


sys.path.append('../../')
#from FenicsCode.Consts import Z0, c0
#from FenicsCode.Utilities.MeshGenerators import get_centred_cube
#import FenicsCode.Utilities.Optimization
#from FenicsCode.PostProcessing import surface_ntff, variational_ntff


# Enable dolfin's form optimizations
sucemfem.Utilities.Optimization.set_dolfin_optimisation()
#FenicsCode.Utilities.Optimization.set_dolfin_optimisation()

# Near-field of an infintesimal dipole
#fname = 'data/dofs-2-0.299792458-0.0166551365556-0.0749481145-1000.pickle'
fname = 'data/dofs-2-0.899377374-0.0499654096667-0.6895226534-1000.pickle'
theta_deg = N.linspace(0, 180, 181)
no_ff_pts = len(theta_deg)
phi_deg = N.zeros(no_ff_pts)
data = pickle.load(open(fname))
lam = c0/data['freq']

k0 = data['freq']*2*N.pi/c0