This is the developmental file for solving the Maxwell's Equation semi-analytically

In [6]:

    

import os
import sys
sys.path.insert(0, os.path.abspath('../../'))

from matplotlib import pyplot as pl
from scipy.integrate import odeint
import numpy as np
import h5py

from dg_maxwell import params

pl.rcParams['figure.figsize'  ] = 9.6, 6.
pl.rcParams['figure.dpi'      ] = 300
pl.rcParams['image.cmap'      ] = 'jet'
pl.rcParams['lines.linewidth' ] = 1.5
pl.rcParams['font.family'     ] = 'serif'
pl.rcParams['font.weight'     ] = 'bold'
pl.rcParams['font.size'       ] = 20
pl.rcParams['font.sans-serif' ] = 'serif'
pl.rcParams['text.usetex'     ] = True
pl.rcParams['axes.linewidth'  ] = 1.5
pl.rcParams['axes.titlesize'  ] = 'medium'
pl.rcParams['axes.labelsize'  ] = 'medium'
pl.rcParams['xtick.major.size'] = 8
pl.rcParams['xtick.minor.size'] = 4
pl.rcParams['xtick.major.pad' ] = 8
pl.rcParams['xtick.minor.pad' ] = 8
pl.rcParams['xtick.color'     ] = 'k'
pl.rcParams['xtick.labelsize' ] = 'medium'
pl.rcParams['xtick.direction' ] = 'in'
pl.rcParams['ytick.major.size'] = 8
pl.rcParams['ytick.minor.size'] = 4
pl.rcParams['ytick.major.pad' ] = 8
pl.rcParams['ytick.minor.pad' ] = 8
pl.rcParams['ytick.color'     ] = 'k'
pl.rcParams['ytick.labelsize' ] = 'medium'
pl.rcParams['ytick.direction' ] = 'in'


    

In [7]:

    

def dydt_0(y, t0):
    '''
    '''
    E_00, B_01 = y
    dydt = [-B_01 * 2 * np.pi, E_00 * 2 * np.pi]
    
    return dydt

def dydt_1(y, t0):
    '''
    '''
    E_01, B_00 = y
    dydt = [B_00 * 2 * np.pi, -E_01 * 2 * np.pi]
    
    return dydt


    

In [8]:

    

E_00_B_01_init = [1.0, 1.]
E_01_B_00_init = [1.0, 1.]

t = np.arange(0., 2., params.delta_t)


    

In [9]:

    

E_00_B_01 = odeint(dydt_0, E_00_B_01_init, t)
E_01_B_00 = odeint(dydt_1, E_01_B_00_init, t)


    

In [10]:

    

print(E_00_B_01)


    

    




[[ 1.          1.        ]
 [ 0.97965143  1.01994268]
 [ 0.95890517  1.03947144]
 ..., 
 [ 1.05353323  0.94343432]
 [ 1.03431333  0.96446697]
 [ 1.01467361  0.98510815]]

In [11]:

    

print(params.delta_t * i * 5)


    

    




---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-11-7ad5fe378132> in <module>()
----> 1 print(params.delta_t * i * 5)

NameError: name 'i' is not defined

In [12]:

    

for i in np.arange(126):
#     i = 20
    print(i)
    E_z = E_00_B_01[5 * i][0] * np.sin(2 * np.pi * params.element_LGL) + E_01_B_00[5 * i][0] * np.cos(2 * np.pi * params.element_LGL)
    B_y = E_01_B_00[5 * i][1] * np.sin(2 * np.pi * params.element_LGL) + E_00_B_01[5 * i][1] * np.cos(2 * np.pi * params.element_LGL)

    f, (ax1, ax2) = pl.subplots(2, sharex = True, sharey = True)

    h5py_data = h5py.File('../../results/hdf5_%02d/dump_timestep_%06d' %(int(params.N_LGL), int(5 * i)) + '.hdf5', 'r')
    E_z_LGL     = (h5py_data['E_z'][:])
    B_y_LGL     = (h5py_data['B_y'][:])

    ax1.plot(params.element_LGL, E_z)
    ax2.plot(params.element_LGL, B_y)

    ax1.plot(params.element_LGL, E_z_LGL, '.')
    ax2.plot(params.element_LGL, B_y_LGL, '.')


    pl.xlabel(r'$x$')
    ax1.set_ylabel(r'$E_z$')
    ax2.set_ylabel(r'$B_y$')
    pl.xlim(-1, 1)
    pl.ylim(-2, 2)
    pl.title(r'$E_z$ Time = %.2f' % (i * 5 * params.delta_t))
    f.savefig('../../results/1D_Wave_images_/%04d' %(i) + '.png')
    pl.close('all')

#     pl.show()


    

    




0
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125






    




---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
<ipython-input-12-0c0ccdd50ca5> in <module>()
      2 #     i = 20
      3     print(i)
----> 4     E_z = E_00_B_01[5 * i][0] * np.sin(2 * np.pi * params.element_LGL) + E_01_B_00[5 * i][0] * np.cos(2 * np.pi * params.element_LGL)
      5     B_y = E_01_B_00[5 * i][1] * np.sin(2 * np.pi * params.element_LGL) + E_00_B_01[5 * i][1] * np.cos(2 * np.pi * params.element_LGL)
      6 

IndexError: index 625 is out of bounds for axis 0 with size 624

In [ ]: