Commit c052c6bee9d8f48a0033b171850a795515d5f698
- Diff rendering mode:
- inline
- side by side
fdtd_bpml_plane.py
(52 / 7)
| 9 | 9 | # program start | |
| 10 | 10 | ||
| 11 | 11 | new_efield = grid.Plane(shape=(301,301), timestep=1, spacestep=2*3*10**8) | |
| 12 | new_efield.eps[:180,0:] *= 2 | ||
| 12 | #new_efield.eps[:180,0:] *= 2 | ||
| 13 | 13 | new_efield.sigmae = grid.PlaneBase(new_efield.shape) | |
| 14 | 14 | new_efield.sigmah = grid.PlaneBase(new_efield.shape) | |
| 15 | ### | ||
| 16 | new_efield.sigmae.y[0,0:301] = new_efield.sigmae.x[300,0:301] = 1000 | ||
| 15 | 17 | ||
| 18 | new_efield.sigmae.y[1,1:300] = new_efield.sigmae.x[299,1:300] = 870 | ||
| 19 | |||
| 20 | new_efield.sigmae.y[2,2:299] = new_efield.sigmae.y[298,2:299] = 660 | ||
| 21 | |||
| 22 | new_efield.sigmae.y[3,3:298] = new_efield.sigmae.y[297,3:298] = 540 | ||
| 23 | |||
| 24 | new_efield.sigmae.y[4,4:297] = new_efield.sigmae.y[296,4:297] = 430 | ||
| 25 | |||
| 26 | new_efield.sigmae.y[5,5:296] = new_efield.sigmae.y[295,5:296] = 350 | ||
| 27 | |||
| 28 | new_efield.sigmae.y[6,6:295] = new_efield.sigmae.y[294,6:295] = 270 | ||
| 29 | |||
| 30 | new_efield.sigmae.y[7,7:294] = new_efield.sigmae.y[293,7:294] = 180 | ||
| 31 | |||
| 32 | new_efield.sigmae.y[8,8:293] = new_efield.sigmae.y[292,8:293] = 120 | ||
| 33 | |||
| 34 | new_efield.sigmae.y[9,9:292] = new_efield.sigmae.y[291,9:292] = 30 | ||
| 35 | |||
| 36 | ### | ||
| 37 | new_efield.sigmae.y[0:301,0] = 5 | ||
| 38 | new_efield.sigmae.y[0:301,300] = 5 | ||
| 39 | |||
| 40 | new_efield.sigmae.y[1:300,1] = 4.5 | ||
| 41 | new_efield.sigmae.y[1:300,299] = 4.5 | ||
| 42 | |||
| 43 | new_efield.sigmae.y[2:299,2] = 4.0 | ||
| 44 | new_efield.sigmae.y[2:299,298] = 4.0 | ||
| 45 | |||
| 46 | new_efield.sigmae.y[3:298,3] = 3.5 | ||
| 47 | new_efield.sigmae.y[3:298,297] = 3.5 | ||
| 48 | |||
| 49 | new_efield.sigmae.y[4:297,4] = 3.0 | ||
| 50 | new_efield.sigmae.y[4:297,296] = 3.0 | ||
| 51 | |||
| 52 | |||
| 53 | ### | ||
| 54 | pylab.imshow(numpy.array(new_efield.sigmae.x + new_efield.sigmae.y, dtype="float")) | ||
| 55 | pylab.savefig("field-sigma-sum.png") | ||
| 56 | pylab.clf() | ||
| 57 | ### | ||
| 58 | |||
| 16 | 59 | new_hfield = grid.Plane(shape=(301,301), timestep=1, spacestep=2*3*10**8) | |
| 17 | new_hfield.eps[:180,0:] *= 2 | ||
| 60 | # new_hfield.eps[:180,0:] *= 2 | ||
| 18 | 61 | new_hfield.zx = numpy.zeros(new_hfield.shape) | |
| 19 | 62 | new_hfield.zy = numpy.zeros(new_hfield.shape) | |
| 20 | new_hfield.sigmae = grid.PlaneBase(new_hfield.shape) | ||
| 21 | new_hfield.sigmah = grid.PlaneBase(new_hfield.shape) | ||
| 63 | # new_hfield.sigmae = grid.PlaneBase(new_hfield.shape) | ||
| 64 | new_hfield.sigmae = new_efield.sigmae | ||
| 65 | # new_hfield.sigmah = grid.PlaneBase(new_hfield.shape) | ||
| 66 | new_hfield.sigmah = new_efield.sigmah | ||
| 22 | 67 | ||
| 23 | 68 | ||
| 24 | 69 | for step in range(1,250): | |
| 25 | 70 | new_efield.y[149:152,149:152] += float(mpmath.sin(step*mpmath.pi/6.125)) | |
| 26 | 71 | new_efield = update_efield( new_efield, new_hfield ) | |
| 27 | new_hfield.zx[150,150] += float(mpmath.sin(step*mpmath.pi/6.125))/2 | ||
| 28 | new_hfield.zy[150,150] += float(mpmath.sin(step*mpmath.pi/6.126))/2 | ||
| 72 | new_hfield.zx[280,150] += float(mpmath.sin(step*mpmath.pi/6.125))/2 | ||
| 73 | new_hfield.zy[280,150] += float(mpmath.sin(step*mpmath.pi/6.126))/2 | ||
| 29 | 74 | new_hfield = update_hfield( new_efield, new_hfield ) | |
| 30 | 75 | ||
| 31 | 76 | h = numpy.array(new_hfield.z, dtype="float") | |
| … | … | ||
| 80 | 80 | if (i+j)%2 == 1: | |
| 81 | 81 | h[i,j] = ( h[i+1,j] + h[i-1,j] + h[i,j+1] + h[i,j-1] ) / 4 | |
| 82 | 82 | ||
| 83 | pylab.imshow( h, interpolation="quadric", norm=matplotlib.colors.Normalize(-1,1.8,True) ) | ||
| 83 | pylab.imshow( h, interpolation="quadric", norm=matplotlib.colors.Normalize(-0.3,0.3,True) ) | ||
| 84 | 84 | pylab.colorbar() | |
| 85 | 85 | ||
| 86 | 86 | pylab.savefig("result/fdtd_plane-%.3d.png" % step) |
