diff --git a/README.md b/README.md index a4329696e..4d3763485 100644 --- a/README.md +++ b/README.md @@ -7,7 +7,7 @@ **Meep** is a free finite-difference time-domain (FDTD) simulation software package to model electromagnetic systems. Meep is an acronym which officially stands for *MIT Electromagnetic Equation Propagation*. Its features include: -- **Free software** under the [GNU GPL](https://en.wikipedia.org/wiki/GNU_General_Public_License). +- **Free and open-source software** under the [GNU GPL](https://en.wikipedia.org/wiki/GNU_General_Public_License). - Complete **scriptability** via [Python](http://meep.readthedocs.io/en/latest/Python_Tutorials/Basics/), [Scheme](http://meep.readthedocs.io/en/latest/Scheme_Tutorials/Basics), or [C++](C++). - Simulation in **1d, 2d, 3d**, and **cylindrical** coordinates. - Distributed memory **parallelism** on any system supporting the [MPI](https://en.wikipedia.org/wiki/MPI) standard. Portable to any Unix-like operating system such as [Linux](https://en.wikipedia.org/wiki/Linux) and [macOS](https://en.wikipedia.org/wiki/macOS). diff --git a/doc/docs/Acknowledgements.md b/doc/docs/Acknowledgements.md index b99d4b371..77454779b 100644 --- a/doc/docs/Acknowledgements.md +++ b/doc/docs/Acknowledgements.md @@ -10,7 +10,7 @@ Meep originated as part of graduate research at MIT with initial contributions b Referencing ----------- -We kindly request that you cite the following reference in any publication for which you use Meep: +We kindly request that you cite the following reference publication in any work for which you use Meep: - A.F. Oskooi, D. Roundy, M. Ibanescu, P. Bermel, J.D. Joannopoulos, and S.G. Johnson, [MEEP: A flexible free-software package for electromagnetic simulations by the FDTD method](http://dx.doi.org/doi:10.1016/j.cpc.2009.11.008), Computer Physics Communications, vol. 181, pp. 687-702 (2010). ([pdf](http://ab-initio.mit.edu/~oskooi/papers/Oskooi10.pdf)) diff --git a/doc/docs/C++_Tutorial.md b/doc/docs/C++_Tutorial.md index 6628cf137..beec08385 100644 --- a/doc/docs/C++_Tutorial.md +++ b/doc/docs/C++_Tutorial.md @@ -2,7 +2,7 @@ # C++ Tutorial --- -Instead of using the [Scheme interface](Scheme_User_Interface/), Meep is also callable as a C++ library by writing a C++ program that links to it. The C++ interface provides the most flexibility in setting up simulations. There are numerous examples in the `tests/` directory of the Meep codebase which cover a wide range of Meep's functionality. These are a good additional reference. Finally, we should also note that, while Meep is nominally in C++, it is perhaps better described as "C+". That is, most of the coding style is C-like with a few C++ features. +Instead of using the [Python interface](Python_User_Interface/), Meep is also callable as a C++ library by writing a C++ program that links to it. The C++ interface provides the most flexibility in setting up simulations. There are numerous examples in the `tests/` directory of the Meep codebase which cover a wide range of Meep's functionality. These are a good additional reference. Finally, we should also note that, while Meep is nominally in C++, it is perhaps better described as "C+". That is, most of the coding style is C-like with a few C++ features. [TOC] diff --git a/doc/docs/Exploiting_Symmetry.md b/doc/docs/Exploiting_Symmetry.md index 414d1a45a..41f53cb4a 100644 --- a/doc/docs/Exploiting_Symmetry.md +++ b/doc/docs/Exploiting_Symmetry.md @@ -8,7 +8,7 @@ An important point to understand is that, when you specify a symmetry, it must b **Meep does not check whether the symmetry is obeyed**. If you specify a symmetry that does not preserve your structure/sources, then the results are undefined. -For the Scheme syntax to specify a symmetry, see the [symmetry reference](Scheme_User_Interface.md#symmetry). There are also examples in [Tutorial/Basics](Python_Tutorials/Basics/#exploiting-symmetry). +For the Python syntax to specify a symmetry, see the [symmetry reference](Python_User_Interface.md#symmetry). There are also examples in [Tutorial/Basics](Python_Tutorials/Basics/#exploiting-symmetry). [TOC] diff --git a/doc/docs/FAQ.md b/doc/docs/FAQ.md index 93492c4c3..f46b7698e 100644 --- a/doc/docs/FAQ.md +++ b/doc/docs/FAQ.md @@ -74,13 +74,13 @@ Usage ### Is there a Python interface? -An official Python interface is currently under development and will likely be released by the end of 2017. An unofficial [Python interface](https://github.com/FilipDominec/python-meep-utils) has been independently developed by researchers at the Institute of Physics at the Czech Academy of Sciences and Ghent University. Unfortunately, this interface has a number of shortcomings including missing support for geometric objects, lack of high-level abstractions for low-level functionality, and limited documentation. These will all be addressed in the official interface. +An official Python interface is currently under development and will be released in early 2018. An unofficial [Python interface](https://github.com/FilipDominec/python-meep-utils) has been developed independently by researchers at the Institute of Physics at the Czech Academy of Sciences and Ghent University. Unfortunately, this interface has a number of shortcomings including missing support for geometric objects, lack of high-level abstractions for low-level functionality, and limited documentation. These will all be addressed in the official interface. ### Why doesn't turning off subpixel averaging work? By default, when Meep assigns a dielectric constant $\varepsilon$ or $\mu$ to each pixel, it uses a carefully designed average of the $\varepsilon$ values within that pixel. This subpixel averaging generally improves the accuracy of the simulation — perhaps counter-intuitively, for geometries with discontinous $\varepsilon$ it is *more* accurate (i.e. closer to the exact Maxwell result for the *discontinuous* case) to do the simulation with the subpixel-averaged (*smoothed*) $\varepsilon$, as long as the averaging is done properly. For details, see the [reference publication](Acknowledgements/#referencing). -Still, there are times when, for whatever reason, you might not want this feature. For example, if your accuracy is limited by other issues, or if you want to skip the wait at the beginning of the simulation for it do to the averaging. In this case, you can disable the subpixel averaging by setting `eps_averaging = False` via the `Simulation` instance (Python) or `(set! eps-averaging? false)` (Scheme) in your script file. See the [User Interface](Python_User_Interface.md). +Still, there are times when, for whatever reason, you might not want this feature. For example, if your accuracy is limited by other issues, or if you want to skip the wait at the beginning of the simulation for it do to the averaging. In this case, you can disable the subpixel averaging by setting `eps_averaging = False` via the `Simulation` class (Python) or `(set! eps-averaging? false)` (Scheme) in your script file. See the [User Interface](Python_User_Interface.md). Even if you disable the subpixel averaging, however, when you output the dielectric function to a file and plot it, you may notice that there are some pixels with intermediate $\varepsilon$ values, right at the boundary between two materials. This has a completely different source. Internally, Meep's simulation is performed on a [Yee grid](Yee_Lattice.md), in which every field component is stored on a slightly different grid which are offset from one another by half-pixels, and the $\varepsilon$ values are also stored on this Yee grid. For output purposes, however, it is more user-friendly to output all fields etcetera on the same grid at the center of each pixel, so all quantities are interpolated onto this grid for output. Therefore, even though the internal $\varepsilon$ values are indeed discontinuous when you disable subpixel averaging, the *output* file will still contain some "averaged" values at interfaces due to the interpolation from the Yee grid to the center-pixel grid. diff --git a/doc/docs/Field_Function_Examples.md b/doc/docs/Field_Function_Examples.md index 3731f5a36..d15866e43 100644 --- a/doc/docs/Field_Function_Examples.md +++ b/doc/docs/Field_Function_Examples.md @@ -2,32 +2,54 @@ # Field Function Examples --- -As described in the [User Interface](Scheme_User_Interface.md), Meep provides several routines to integrate, analyze, and output arbitrary user-specified functions of the field components. See the functions whose names end with `-field-function`. This facility, while powerful, requires a bit more Scheme programming than most Meep usage, and is best illustrated by a few examples. +As described in the [User Interface](Python_User_Interface.md), Meep provides several routines to integrate, analyze, and output arbitrary user-specified functions of the field components. See the functions whose names end with `_field_function`. This facility, while powerful, requires a bit more programming than most Meep usage, and is best illustrated by a few examples. Every field-function that can be passed to these routines is of the form *f*(**r**,components...), where **r** is a position vector and "components..." are zero or more field components that the function depends on. The set of desired components is user-specified. As an arbitrary example, suppose we are interested in the strange function: $$f(\mathbf{r}, E_x, H_z, \varepsilon) = x |\mathbf{r}| + E_x - \varepsilon H_z$$ -We would define this function, in Scheme, by: +We would define this function by: + +*Python* +```py +def f(r, ex, hz, eps): + return (r.x * r.norm() + ex) - (eps * hz) +``` + +*Scheme* ```scm (define (f r ex hz eps)    (- (+ (* (vector3-x r) (vector3-norm r)) ex) (* eps hz))) ``` -Note that the `r` argument is a `vector3`, and can be manipulated by the functions defined in the [Libctl User Reference](https://libctl.readthedocs.io/en/latest/Libctl_User_Reference/). +Note that the `r` argument is a `Vector3` (Python) or `vector3` (Scheme), and can be manipulated by the functions defined in the [Libctl User Reference](https://libctl.readthedocs.io/en/latest/Libctl_User_Reference/). -Now, suppose we want to compute the integral of this function, over the whole computational cell. We can do this by calling the function `integrate-field-function`, as follows: +Now, suppose we want to compute the integral of this function, over the whole computational cell. We can do this by calling the function `integrate_field_function` (Python) or `integrate-field-function` (Scheme), as follows: +*Python* +```py +print("The integral of our weird function is: {}" + .format(meep.Simulation.integrate_field_function([meep.Ex, meep.Hz, meep.Dielectric], f))) +``` + +*Scheme* ```scm -(print "The integral of our weird function is: " -       (integrate-field-function (list Ex Hz Dielectric) f) "\n") +(print "The integral of our weird function is: " + (integrate-field-function (list Ex Hz Dielectric) f) "\n") ``` -Note that the first argument to `integrate-field-function` is a `list` (a standard Scheme type) of `component` constants, specifying in order the list of field components our function `f` expects to be passed. Meep will then call `f` for every point in the computational cell in parallel on a parallel machine, and return the integral approximated by a [trapezoidal rule](https://en.wikipedia.org/wiki/trapezoidal_rule). +Note that the first argument to `integrate_field_function` (Python) or `integrate-field-function` (Scheme) is a list (a standard Python/Scheme type) of `component` constants, specifying in order the list of field components our function `f` expects to be passed. Meep will then call `f` for every point in the computational cell in parallel on a parallel machine, and return the integral approximated by a [trapezoidal rule](https://en.wikipedia.org/wiki/trapezoidal_rule). -You can also specify an optional third argument to `integrate-field-function`, specifying an integration volume in case you don't want the integral over the whole computational cell. For example, the following code computes the integral of `f` along a line from (-1,0,0) to (1,0,0): +You can also specify an optional third argument to `integrate_field_function` or `integrate-field-function`, specifying an integration volume in case you don't want the integral over the whole computational cell. For example, the following code computes the integral of `f` along a line from (-1,0,0) to (1,0,0): +*Python* +```py +print("The integral of our weird function from (-1,0,0) to (1,0,0) is: {}" + .format(meep.Simulation.integrate_field_function([meep.Ex, meep.Hz, meep.Dielectric], f, meep.Volume(size=meep.Vector3(1), center=meep.Vector3())))) +``` + +*Scheme* ```scm (print "The integral of our weird function from (-1,0,0) to (1,0,0) is: "        (integrate-field-function (list Ex Hz Dielectric) f (volume (size 1 0 0) (center 0 0 0))) "\n") @@ -35,6 +57,13 @@ You can also specify an optional third argument to `integrate-field-function`, s Instead of computing the integral, Meep also provides a function to compute the maximum absolute value of our given function: +*Python* +```py +print("The maximum absolute value of our weird function from (-1,0,0) to (1,0,0) is: {}" + .format(meep.Simulation.max_abs_field_function([meep.Ex, meep.Hz, meep.Dielectric], f, meep.Volume(size=meep.Vector3(1), center=meep.Vector3())))) +``` + +*Scheme* ```scm (print "The maximum absolute value of our weird function from (-1,0,0) to (1,0,0) is: "        (max-abs-field-function (list Ex Hz Dielectric) f (volume (size 1 0 0) (center 0 0 0))) "\n") @@ -42,13 +71,19 @@ Instead of computing the integral, Meep also provides a function to compute the Finally, we can also output our function to an HDF5 file, similar to the built-in functions to output selected field components, and so on. The following outputs an HDF5 file consisting of our function `f` evaluated at every point in the computational cell: +*Python* +```py +meep.Simulation.output_field_function("weird-function", [meep.Ex, meep.Hz, meep.Dielectric], f) +``` + +*Scheme* ```scm (output-field-function "weird-function" (list Ex Hz Dielectric) f) ``` -Here, the first argument is used for the name of the dataset within the HDF5, and is also used for the name of the HDF5 file itself plus a `.h5` suffix and a time stamp, unless you have specified the output file via `to-appended` or other means. +Here, the first argument is used for the name of the dataset within the HDF5, and is also used for the name of the HDF5 file itself plus a `.h5` suffix and a time stamp, unless you have specified the output file via `to_appended` or `to-appended` or other means. -The above example calls the integration, maximum, and output routines only once, at the current time. Often, you will want to pass them to `run-until` instead, using `at-every` to print or output at periodic time intervals. A common mistake is to do something like the following: +The above example calls the integration, maximum, and output routines only once, at the current time. Often, you will want to pass them to `meep.Simulation.run(..., until=...)` (Python) or `run-until` (Scheme) instead, using `at_every` or `at-every` to print or output at periodic time intervals. In Scheme, a common mistake is to do something like the following: ```scm (run-until 200 (at-every 1 (output-field-function "weird-function" (list Ex Hz Dielectric) f))) @@ -69,4 +104,4 @@ As described in [Synchronizing the Magnetic and Electric Fields](Synchronizing_t (run-until 200 (synchronized-magnetic (at-every 1 my-weird-output))) ``` -See also the section [writing your own step functions](Scheme_User_Interface.md#writing-your-own-step-functions). +See also the section "Writing Your Own Step Functions" in the [Python](Python_User_Interface.md#writing-your-own-step-functions) or [Scheme](Scheme_User_Interface.md#writing-your-own-step-functions) interface. diff --git a/doc/docs/License_and_Copyright.md b/doc/docs/License_and_Copyright.md index 0c3685798..ed5f8c858 100644 --- a/doc/docs/License_and_Copyright.md +++ b/doc/docs/License_and_Copyright.md @@ -2,7 +2,7 @@ # License and Copyright --- -Meep is copyright © 2005–2017, Massachusetts Institute of Technology. +Meep is copyright © 2005–2018, Massachusetts Institute of Technology. Meep is free software. You can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2 of the License, or at your option any later version. diff --git a/doc/docs/Materials.md b/doc/docs/Materials.md index 9bc3f1644..d674918d5 100644 --- a/doc/docs/Materials.md +++ b/doc/docs/Materials.md @@ -8,7 +8,7 @@ However, $\varepsilon$ is not only a function of position. In general, it also d Similarly for the relative permeability $\mu$(**r**), for which dispersion, nonlinearity, and anisotropy are all supported as well. -In this section, we describe the form of the equations and material properties that Meep can simulate. The actual user interface where these properties are specified is described in the [User Interface](Scheme_User_Interface.md). +In this section, we describe the form of the equations and material properties that Meep can simulate. The actual user interface where these properties are specified is described in the [User Interface](Python_User_Interface.md). [TOC] @@ -80,7 +80,7 @@ Often, you only care about the absorption loss in a narrow bandwidth, where you One approach to this problem would be allowing you to specify a constant, frequency-independent, imaginary part of $\varepsilon$, but this has the disadvantage of requiring the simulation to employ complex fields which double the memory and time requirements, and also tends to be numerically unstable. Instead, the approach in Meep is for you to set the conductivity $\sigma_D$ (or $\sigma_B$ for an imaginary part of $\mu$), chosen so that $\mathrm{Im}\, \varepsilon = \varepsilon_\infty \sigma_D / \omega$ is the correct value at your frequency $\omega$ of interest. Note that, in Meep, you specify $f = \omega/2\pi$ instead of $\mu$ for the frequency, however, so you need to include the factor of 2$\pi$ when computing the corresponding imaginary part of $\varepsilon$. Conductivities can be implemented with purely real fields, so they are not nearly as expensive as implementing a frequency-independent complex $\varepsilon$ or $\mu$. -For example, suppose you want to simulate a medium with $\varepsilon = 3.4 + 0.101i$ at a frequency 0.42 (in your Meep units), and you only care about the material in a narrow bandwidth around this frequency (i.e. you don't need to simulate the full experimental frequency-dependent permittivity). Then, in Meep, you could use `(make medium (epsilon 3.4) (D-conductivity (/ (* 2 pi 0.42 0.101) 3.4))`; i.e. $\varepsilon_\infty = \mathrm{Re}\,\varepsilon = 3.4$ and $\sigma_D = \omega \mathrm{Im} \varepsilon / \varepsilon_\infty = (2\pi 0.42) 0.101 / 3.4$. +For example, suppose you want to simulate a medium with $\varepsilon = 3.4 + 0.101i$ at a frequency 0.42 (in your Meep units), and you only care about the material in a narrow bandwidth around this frequency (i.e. you don't need to simulate the full experimental frequency-dependent permittivity). Then, in Meep, you could use `meep.Medium(epsilon=3.4, D_conductivity=2*math.pi*0.42*0.101/3.4)`; i.e. $\varepsilon_\infty = \mathrm{Re}\,\varepsilon = 3.4$ and $\sigma_D = \omega \mathrm{Im} \varepsilon / \varepsilon_\infty = (2\pi 0.42) 0.101 / 3.4$. **Note**: the "conductivity" in Meep is slightly different from the conductivity you might find in a textbook, because for computational convenience it appears as $\sigma_D \mathbf{D}$ in our Maxwell equations rather than the more-conventional $\sigma \mathbf{E}$; this just means that our definition is different from the usual electric conductivity by a factor of $\varepsilon$. Also, just as Meep uses the dimensionless relative permittivity for $\varepsilon$, it uses nondimensionalized units of 1/*a* (where *a* is your unit of distance) for the conductivities $\sigma_{D,B}$. If you have the electric conductivity $\sigma$ in SI units and want to convert to $\sigma_D$ in Meep units, you can simply use the formula: $\sigma_D = (a/c) \sigma / \varepsilon_r \varepsilon_0$ where *a* is your unit of distance in meters, *c* is the vacuum speed of light in m/s, $\varepsilon_0$ is the SI vacuum permittivity, and $\varepsilon_r$ is the real relative permittivity. diff --git a/doc/docs/Parallel_Meep.md b/doc/docs/Parallel_Meep.md index 01ffb415c..190c2307f 100644 --- a/doc/docs/Parallel_Meep.md +++ b/doc/docs/Parallel_Meep.md @@ -22,11 +22,17 @@ The parallel version of Meep is designed to operate completely transparently: yo In order to run MPI programs, you typically have to use a command like `mpirun` with an argument to say how many processes you want to use. Consult your MPI documentation. For example, with many popular MPI implementations, to run with 4 processes you would use something like: +*Python* +```sh +mpirun -np 4 python foo.py > foo.out +``` + +*Scheme* ```sh mpirun -np 4 meep-mpi foo.ctl > foo.out ``` -There is one important requirement: every MPI process must be able to read the `foo.ctl` input file or whatever your control file is called. On most systems, this is no problem, but if for some reason your MPI processes don't all have access to the local filesystem then you may need to make copies of your input file or something. +There is one important requirement: every MPI process must be able to read the `foo.py` or `foo.ctl` input file or whatever your control file is called. On most systems, this is no problem, but if for some reason your MPI processes don't all have access to the local filesystem then you may need to make copies of your input file or something. You cannot run Meep interactively on multiple processors. @@ -39,15 +45,15 @@ However, there is an alternative strategy for parallelization. If you have many Technical Details ----------------- -When you run Meep under MPI, the following is a brief description of what is happening behind the scenes. For the most part, you shouldn't *need* to know this stuff. Just use the same `.ctl` exactly as you would for a uniprocessor simulation. +When you run Meep under MPI, the following is a brief description of what is happening behind the scenes. For the most part, you shouldn't *need* to know this stuff. Just use the same `.py` or `.ctl` file exactly as you would for a uniprocessor simulation. -First, every MPI process executes the `.ctl` file in parallel. The processes communicate however, to only perform one simulation in sync with one another. In particular, the computational cell is divided into "chunks", one per process, to roughly equally divide the work and the memory. +First, every MPI process executes the `.py` or `.ctl` file in parallel. The processes communicate however, to only perform one simulation in sync with one another. In particular, the computational cell is divided into "chunks", one per process, to roughly equally divide the work and the memory. -When you time-step (via `run-until` or whatever), the chunks are time-stepped in parallel, communicating the values of the pixels on their boundaries with one another. In general, any Meep function that performs some collective operation over the whole computational cell or a large portion thereof is parallelized, including: time-stepping, HDF5 I/O, accumulation of flux spectra, and field integration via `integrate-field-function`, although the *results* are communicated to all processes. +When you time-step [via `meep.Simulation.run(until=...)` or `run-until` or whatever], the chunks are time-stepped in parallel, communicating the values of the pixels on their boundaries with one another. In general, any Meep function that performs some collective operation over the whole computational cell or a large portion thereof is parallelized, including: time-stepping, HDF5 I/O, accumulation of flux spectra, and field integration via `integrate_field_function` (Python) or `integrate-field-function` (Scheme), although the *results* are communicated to all processes. -Computations that only involve isolated points, such as `get-field-point` or `harminv` analyses, are performed by all processes redundantly. In the case of `get-field-point`, Meep figures out which process "knows" the field at the given field, and then sends the field value from that process to all other processes. This is harmless because such computations are rarely performance bottlenecks. +Computations that only involve isolated points, such as `get_field_point` (Python) or `get-field-point` (Scheme), or `Harminv` (Python) or `harminv` (Scheme) analyses, are performed by all processes redundantly. In the case of `get_field_point` or `get-field-point`, Meep figures out which process "knows" the field at the given field, and then sends the field value from that process to all other processes. This is harmless because such computations are rarely performance bottlenecks. -Although all processes execute the `.ctl` in parallel, `print` statements are ignored for all process but one (process \#0). In this way, you only get one copy of the output. +Although all processes execute the `.py` or `.ctl` file in parallel, `print` statements are ignored for all process but one (process \#0). In this way, you only get one copy of the output. If for some reason you need to distinguish different MPI processes in your `.ctl` file, you can use the following two functions: diff --git a/doc/docs/Scheme_Tutorials/Basics.md b/doc/docs/Scheme_Tutorials/Basics.md index 21a3394cb..5c6010790 100644 --- a/doc/docs/Scheme_Tutorials/Basics.md +++ b/doc/docs/Scheme_Tutorials/Basics.md @@ -17,7 +17,7 @@ Don't worry, though — simple things are simple and you don't need to be an The ctl file is actually implemented on top of the [libctl](https://libctl.readthedocs.io) library, a set of utilities that are in turn built on top of the Scheme language. Thus, there are three sources of possible commands and syntax for a ctl file: -- [Scheme](https://en.wikipedia.org/wiki/Scheme_programming_language) a powerful and beautiful programming language developed at MIT. The syntax is particularly simple: all statements are of the form `(function` `arguments...)`. We run Scheme under the [Guile](https://en.wikipedia.org/wiki/GNU_Guile) interpreter which is designed to be plugged into programs as a scripting and extension language. You don't need to know much Scheme for a basic ctl file, but it is always there if you need it. More information is available on [Guile and Scheme](Guile_and_Scheme_Information.md). +- [Scheme](https://en.wikipedia.org/wiki/Scheme_programming_language) a powerful and beautiful programming language developed at MIT. The syntax is particularly simple: all statements are of the form `(function` `arguments...)`. We run Scheme under the [Guile](https://en.wikipedia.org/wiki/GNU_Guile) interpreter which is designed to be plugged into programs as a scripting and extension language. You don't need to know much Scheme for a basic ctl file, but it is always there if you need it. More information is available on [Guile and Scheme](../Guile_and_Scheme_Information.md). - [libctl](https://libctl.readthedocs.io/), a library that we built on top of Guile to simplify communication between Scheme and scientific computation software. libctl sets the basic tone of the interface and defines a number of useful functions (such as multi-variable optimization, numeric integration, and so on). See the [libctl manual](https://libctl.readthedocs.io). - Meep itself, which defines all the interface features that are specific to FDTD calculations. This manual is primarily focused on documenting these features. @@ -70,13 +70,13 @@ Now that we have the structure, we need to specify the current sources, which is Here, we gave the source a frequency of 0.15, and specified a `continuous-src` which is just a fixed-frequency sinusoid $\exp(-i\omega t)$ that (by default) is turned on at $t=0$. Recall that, in [Meep units](../Introduction.md#units-in-meep), frequency is specified in units of $2\pi c$, which is equivalent to the inverse of vacuum wavelength. Thus, 0.15 corresponds to a vacuum wavelength of about $1/0.15=6.67$, or a wavelength of about 2 in the $\varepsilon=12$ material—thus, our waveguide is half a wavelength wide, which should hopefully make it single-mode. In fact, the cutoff for single-mode behavior in this waveguide is analytically solvable, and corresponds to a frequency of 1/2√11 or roughly 0.15076. Note also that to specify a $J_z$, we specify a component $Ez$ (e.g. if we wanted a magnetic current, we would specify `Hx`, `Hy`, or `Hz`). The current is located at $(-7,0)$, which is 1 unit to the right of the left edge of the cell — we always want to leave a little space between sources and the cell boundaries, to keep the boundary conditions from interfering with them. -Speaking of boundary conditions, we want to add absorbing boundaries around our cell. Absorbing boundaries in Meep are handled by [perfectly matched layers](Perfectly_Matched_Layer.md) (PML) — which aren't really a boundary condition at all, but rather a fictitious absorbing material added around the edges of the cell. To add an absorbing layer of thickness 1 around all sides of the cell, we do: +Speaking of boundary conditions, we want to add absorbing boundaries around our cell. Absorbing boundaries in Meep are handled by [perfectly matched layers](../Perfectly_Matched_Layer.md) (PML) — which aren't really a boundary condition at all, but rather a fictitious absorbing material added around the edges of the cell. To add an absorbing layer of thickness 1 around all sides of the cell, we do: ```scm (set! pml-layers (list (make pml (thickness 1.0)))) ``` -`pml-layers` is a list of `pml` objects — you may have more than one `pml` object if you want PML layers only on certain sides of the cell, e.g. `(make pml (thickness 1.0) (direction X) (side High))` specifies a PML layer on only the $+x$ side. Now, we note an important point: **the PML layer is *inside* the cell**, overlapping whatever objects you have there. So, in this case our PML overlaps our waveguide, which is what we want so that it will properly absorb waveguide modes. The finite thickness of the PML is important to reduce numerical reflections; see [Perfectly Matched Layer](Perfectly_Matched_Layer.md) for more information. +`pml-layers` is a list of `pml` objects — you may have more than one `pml` object if you want PML layers only on certain sides of the cell, e.g. `(make pml (thickness 1.0) (direction X) (side High))` specifies a PML layer on only the $+x$ side. Now, we note an important point: **the PML layer is *inside* the cell**, overlapping whatever objects you have there. So, in this case our PML overlaps our waveguide, which is what we want so that it will properly absorb waveguide modes. The finite thickness of the PML is important to reduce numerical reflections; see [Perfectly Matched Layer](../Perfectly_Matched_Layer.md) for more information. Meep will discretize this structure in space and time, and that is specified by a single variable, `resolution`, that gives the number of pixels per distance unit. We'll set this resolution to 10, which corresponds to around 67 pixels/wavelength, or around 20 pixels/wavelength in the high-dielectric material. In general, at least 8 pixels/wavelength in the highest dielectric is a good idea. This will give us a $160\times80$ cell. @@ -401,7 +401,7 @@ Again, we must run *both* simulations in order to get the normalization right. T Modes of a Ring Resonator ------------------------- -As described in the [Introduction](Introduction.md#resonant-modes), another common task for FDTD simulation is to find the resonant modes — frequencies and decay rates — of some electromagnetic cavity structure. You might want to read that introduction again to recall the basic computational strategy. Here, we will show how this works for perhaps the simplest example of a dielectric cavity: a **ring resonator**, which is simply a waveguide bent into a circle. This can be also found in the `examples/ring.ctl` file included with Meep. In fact, since this structure has cylindrical symmetry, we can simulate it *much* more efficiently [by using cylindrical coordinates](Ring_Resonator_in_Cylindrical_Coordinates.md), but for illustration here we'll just use an ordinary 2d simulation. +As described in the [Introduction](../Introduction.md#resonant-modes), another common task for FDTD simulation is to find the resonant modes — frequencies and decay rates — of some electromagnetic cavity structure. You might want to read that introduction again to recall the basic computational strategy. Here, we will show how this works for perhaps the simplest example of a dielectric cavity: a **ring resonator**, which is simply a waveguide bent into a circle. This can be also found in the `examples/ring.ctl` file included with Meep. In fact, since this structure has cylindrical symmetry, we can simulate it *much* more efficiently [by using cylindrical coordinates](Ring_Resonator_in_Cylindrical_Coordinates.md), but for illustration here we'll just use an ordinary 2d simulation. As before, we'll define some parameters to describe the geometry, so that we can easily change the structure: @@ -508,7 +508,7 @@ which differs by about 0.000001 ($10^{-6}$) from the earlier estimate; the diffe ### Exploiting Symmetry -In this case, because we have a mirror symmetry plane (the $x$ axis) that preserves *both* the structure *and* the sources, we can **exploit this mirror symmetry to speed up the computation**. See also [Exploiting Symmetry](Exploiting_Symmetry.md). In particular, everything about the input file is the same except that we add a single line, right after we specify the `sources`: +In this case, because we have a mirror symmetry plane (the $x$ axis) that preserves *both* the structure *and* the sources, we can **exploit this mirror symmetry to speed up the computation**. See also [Exploiting Symmetry](../Exploiting_Symmetry.md). In particular, everything about the input file is the same except that we add a single line, right after we specify the `sources`: ```scm (set! symmetries (list (make mirror-sym (direction Y)))) diff --git a/doc/docs/Synchronizing_the_Magnetic_and_Electric_Fields.md b/doc/docs/Synchronizing_the_Magnetic_and_Electric_Fields.md index 7060c56a5..0f12c21fa 100644 --- a/doc/docs/Synchronizing_the_Magnetic_and_Electric_Fields.md +++ b/doc/docs/Synchronizing_the_Magnetic_and_Electric_Fields.md @@ -15,18 +15,30 @@ All of this process is handled for you in Meep by a single step function: `synch For example, if you do: +*Python* +```py +meep.Simulation.run(meep.output_poynting,meep.output_tot_pwr,until=200) +``` + +*Scheme* ```scm (run-until 200 output-poynting output-tot-pwr) ``` it outputs the Poynting vector and the total energy density in the electric and magnetic fields at each timestep, but it only does so to first-order accuracy because those computations combine unsynchronized electric and magnetic fields. Instead, if you do +*Python* +```py +meep.Simulation.run(meep.synchronized_magnetic(meep.output_poynting,meep.output_tot_pwr,until=200)) +``` + +*Scheme* ```scm (run-until 200 (synchronized-magnetic output-poynting output-tot-pwr)) ``` it will output the same quantities, but more accurately because the fields will be synchronized. Of course, **there is a price**: synchronizing the fields takes time, and also increases the memory usage in order to backup the unsynchronized fields. -Alternatively, if you want to synchronize the magnetic and electric fields in some context other than that of a step function, e.g. you are doing some computation like `integrate-field-function` outside of the timestepping, you can instead call two lower-level functions. Before doing your computations, you should call `(meep-fields-synchronize-magnetic-fields fields)` to synchronize the magnetic fields with the electric fields, and after your computation you should call `(meep-fields-restore-magnetic-fields fields)` to restore the fields to their unsynchronized state for timestepping. In the C++ interface, these correspond to `fields::synchronize_magnetic_fields` and `fields::restore_magnetic_fields`. If you *don't* call `meep-fields-restore-magnetic-fields` before timestepping, then the fields will be re-synchronized after *every* timestep, which will greatly increase the cost of timestepping. +Alternatively, if you want to synchronize the magnetic and electric fields in some context other than that of a step function, e.g. you are doing some computation like `integrate_field_function` (Python) or `integrate-field-function` (Scheme) outside of the timestepping, you can instead call two lower-level functions. Before doing your computations, you should call `meep.Simulation.fields.synchronize_magnetic_fields()` (Python) or `(meep-fields-synchronize-magnetic-fields fields)` (Scheme) to synchronize the magnetic fields with the electric fields, and after your computation you should call `meep.Simulation.fields.restore_magnetic_fields()` (Python) or `(meep-fields-restore-magnetic-fields fields)` (Scheme) to restore the fields to their unsynchronized state for timestepping. In the C++ interface, these correspond to `fields::synchronize_magnetic_fields` and `fields::restore_magnetic_fields`. If you *don't* call `meep.Simulation.fields.restore_magnetic_fields` or `meep-fields-restore-magnetic-fields` before timestepping, then the fields will be re-synchronized after *every* timestep, which will greatly increase the cost of timestepping. -In future versions, we may decide to synchronize the fields automatically whenever you output something like the Poynting vector or do another field computation that involves both magnetic and electric fields, but currently you must do this manually. In any case, Meep does no additional work when you nest synchronization calls, so it is harmless to insert redundant field synchronizations. The `flux-in-box` and `field-energy-in-box` routines are already automatically synchronized, however. \ No newline at end of file +In future versions, we may decide to synchronize the fields automatically whenever you output something like the Poynting vector or do another field computation that involves both magnetic and electric fields, but currently you must do this manually. In any case, Meep does no additional work when you nest synchronization calls, so it is harmless to insert redundant field synchronizations. The `flux_in_box` (Python) or `flux-in-box` (Scheme) and `field_energy_in_box` (Python) or `field-energy-in-box` (Scheme) routines are already automatically synchronized, however. \ No newline at end of file diff --git a/doc/docs/Units_and_Nonlinearity.md b/doc/docs/Units_and_Nonlinearity.md index 0c6ba1fef..b5d317f64 100644 --- a/doc/docs/Units_and_Nonlinearity.md +++ b/doc/docs/Units_and_Nonlinearity.md @@ -31,6 +31,6 @@ The key to using correct magnitudes is nonlinear calculations in Meep is to real For example, suppose now that we are given some experimental value of $n_2$ in "real" units, say $3\times10^{–8}$μm2/W (silica glass). Of course, this is very small (semiconductors can be much more nonlinear), so to compensate suppose let's plan to use an unrealistic 1 MW of power in our structure (say a waveguide). To implement this in Meep, for convenience we'll use μm as our unit of distance and W as our units of power. First, we'll set $\chi^{(3)}$ from $n_2$ and *n* in these units: $\chi^{(3)} = 4n^2 (3\times 10^{-8})/3$ where *n* is the linear index. Then we simply monitor the power going through our waveguide (or whatever) and change the current amplitude (note that power ∼ $J^2$) until we get $10^6$. -To monitor the power in a structure, we can use a variety of functions. See the [User Interface](Scheme_User_Interface.md). One can get the power flux directly through the `meep-fields-flux-in-box` function. Or, one can alternatively get the intensity at a single point by calling `get-field-point` and passing `Sx` etc. for the component. One thing to be cautious about is that these return the power or intensity at one instant in time, and not the time-average unless you use complex-valued fields (which are problematic for nonlinear systems). +To monitor the power in a structure, we can use a variety of functions. See the [User Interface](Python_User_Interface.md). One can get the power flux directly through the `flux_in_box` function. Or, one can alternatively get the intensity at a single point by calling `get_field_point` and passing `Sx` etc. for the component. One thing to be cautious about is that these return the power or intensity at one instant in time, and not the time-average unless you use complex-valued fields (which are problematic for nonlinear systems). You may have been hoping for a simple formula: set the current to *x* to get *y* power. However, this is not feasible since the amount of power or field intensity you get from a current source depends on the source geometry, the dielectric structure, and so on. And a formula for the units of current is not terribly useful because usually the current source in an FDTD calculation is artifically inserted to create the field, and doesn't correspond to the current source in a physical experiment. diff --git a/doc/docs/index.md b/doc/docs/index.md index 6acc43f6b..c3c1a6989 100644 --- a/doc/docs/index.md +++ b/doc/docs/index.md @@ -42,7 +42,7 @@ Documentation See the navigation sidebar at left. In particular, the [Introduction](Introduction.md) and [Tutorial/Basics](Python_Tutorials/Basics.md) are the most important things to read. There is also an [FAQ](FAQ.md). -Please [cite Meep](Acknowledgements.md#referencing) in any publication for which you found it useful. +Please [cite the reference publication](Acknowledgements.md#referencing) in any work for which you found Meep useful. ### Mailing Lists