Main Lemma Repository
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

lemma.h 14KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331
  1. /* This file is part of Lemma, a geophysical modelling and inversion API */
  2. /* This Source Code Form is subject to the terms of the Mozilla Public
  3. * License, v. 2.0. If a copy of the MPL was not distributed with this
  4. * file, You can obtain one at http://mozilla.org/MPL/2.0/. */
  5. /**
  6. @file
  7. @author Trevor Irons
  8. @date 06/19/2009 09:12:20 AM The Birth of Lemma!
  9. @version $Id: lemma.h 203 2015-01-09 21:19:04Z tirons $
  10. **/
  11. // \image html lemma.png
  12. /** \mainpage Lemma is an ElectroMagnetics Modelling API
  13. \authors Trevor Irons and M. Andrew Kass and others
  14. Originally Lemma was intended as a recursive acronym standing for
  15. <B>L</B>emma is an <B>E</B>lectro<B>M</B>agnetics <B>M</B>odelling <B>A</B>PI.
  16. As the breadth of the project has expanded, the name has remained
  17. appropriate in a more literal sense. Lemma is a flexible cross-platform
  18. library delivering an expressive API that can be used to easily create
  19. versatile programs. Lemma is not itself a program, instead it is a
  20. collection of building blocks to make geophysical applications.
  21. We retain this name because:
  22. - In mathematics a Lemma is a proven proposition which is used as a
  23. stepping stone to a larger result rather than as a statement in-and-of
  24. itself.
  25. - In addition to the electromagnetic modelling, some other facilities are
  26. provided such as numerical optimization and inversion capabilities. These
  27. tools are also considered stepping stones to final products.
  28. We feel that this is a particularly approprate name, as Lemma's
  29. API can be leveraged create powerful applications such as forward
  30. modelling and inverting frequency and time-domain
  31. surveys of arbitrary survey design, sNMR surveys, CSAMT and more.
  32. \section Motivation
  33. Why another Geophysical EM project? For starters, there aren't that many
  34. quality open source packages out there. Those that do exist are generally
  35. specialized to perform a single task and extending them is a major undertaking.
  36. Lemma's approach is much different, by providing a set of general tools users
  37. can easily assemble applications that suite their needs. Furthermore, most are
  38. written in either Fortran or MATLAB, and can be difficult to integrate into
  39. multiphysics applications. To our knowlege, Lemma is the only C++ EM simulation
  40. package freely available.
  41. \section Capabilities Capabilities
  42. In the long term, we have many goals for this software project. Due to its
  43. design, Lemma can be built upon and extended easily. The initial aim is to
  44. provide flexible 1D and 3D EM modelling in the time and frequency domains.
  45. \subsection FDM Frequency-domain forward modelling
  46. Lemma was initially called EMMODFD: Electromagnetic Modelling in the Frequency
  47. Domain. As such this is the most mature area of Lemma.
  48. \par 1D
  49. Frequency domain solutions to electrical and magnetic dipoles can be computed
  50. quasi-analytically in 1D. Calculations can be made in or above the layered
  51. media, and complex electrical conductivity and magnetic susceptibility are
  52. supported according to the Cole-Cole model. Sources may be embdedded in the
  53. media or in the resisitive air layer. Lemma can also can compute fields due
  54. to arbitrarily shaped ungrounded wire loops, topography of the loops is also
  55. supported. Two separate approaches to solving the Hankel transform, one
  56. based on Anderson's digitial filtering technique, and another based on Gaussian
  57. quadrature.
  58. \par 3D
  59. A fast 3D solver that can modify the 1D results based
  60. on arbitrary electrical conductivity model is nearing completion and is
  61. provided in a separate module.
  62. \par future work
  63. We are also planning on supporting grounded wires in the near future.
  64. \subsection TDM Time-domain forward modelling
  65. A 1D time-domain solution has been implemented that utilises both a
  66. dipole source as well as a wire loop. Currently, only one receiver is
  67. modelled at a time, but will be generalised. In addition, utilities
  68. to read in data files for modelling have been implemented.
  69. \subsection DataFormats Data Formats
  70. The EM community is plagued with myriad data formats. Often each equiptment
  71. manufacturer provides their own data format and interoperability is a
  72. constant struggle. We are working on a flexible data format based on the XML
  73. format that can be adapted to many types of data. The template for this
  74. format will be publically released and we hope it catches on in the community.
  75. At the least, it will provide a mechanism to compare datasets and datatypes
  76. within Lemma.
  77. \section Modules Modules
  78. Due to Lemma's design, it is easy to extend the platform. In some cases this
  79. extension results in adding functionality that is not directly related to
  80. ElectroMagnetics. The following modules utilise parts of Lemma to provide
  81. their functionality.
  82. \section Tutorials
  83. - \ref Tutorial - Basic intruduction to Lemma, including aquiring and
  84. compiling the code, class structure, and building your own
  85. applications.
  86. - \ref Extending Tutorial on how to extend Lemma.
  87. \section Development Development and design
  88. Ths package was initially developed by the Center for Gravity, Electrical, and
  89. Magnetic Studies (CGEM) at the Colorado School of Mines (CSM), the United
  90. States Geological Survey (USGS), and Broken Spoke Development, LLC. Where it drew
  91. on work by many others including Ki Ha Lee, and Walt Anderson. All new work and
  92. interfaces are written entirely in C++. Several small external projects are
  93. included, which are written in standard C, and FORTRAN 77. We adapt a
  94. modern, test driven, object oriented, C++ framework.
  95. More recent development has been undertaken at the University of Utah through the Energy
  96. and Geoscience Institute.
  97. \section Legalities
  98. \subsection Copyrights
  99. The following copyrights apply to the source.
  100. Most of the code was developed either by Trevor
  101. Copyright (C) 2008-2010 Trevor Irons <trevor.irons@lemmasoftware.org> or
  102. M. Andrew Kass Copyright (C) 2010 <mkass@usgs.gov>.
  103. The 1D EM solver was derived (but updated heavily) from a Fortran
  104. programme written by Ki Ha Lee in 1984. We have communicated with Ki Ha,
  105. and he assured us that this code is in the public domain.
  106. A Gaussian quadrature hankel transform originally written by Alan Chave was
  107. ported to C++. This code is in the public domain, and the source code was
  108. published in Geophysics.
  109. A digital filtering approach to the Hankel transform written by Walt
  110. Anderson was also rewritten for Lemma. The original Fortran code is also in
  111. the public domain.
  112. Please note that Ki Ha Lee and Walt Anderson had no part in this work, and
  113. the above should not be interpreted as any sort of endorsement by those
  114. parties.
  115. \subsection License
  116. This Source Code Form is subject to the terms of the Mozilla Public
  117. License, v. 2.0. If a copy of the MPL was not distributed with this
  118. file, You can obtain one at http://mozilla.org/MPL/2.0/.
  119. \section Contributing Suggestions and contributions
  120. We welcome contributions and suggestions. Feel free to email the development
  121. team at info@lemmasoftware.org.
  122. Under the terms of the MPL, if you modify a Lemma file, you are obligated to
  123. share those contributions back with the community.
  124. \section Useful Useful links
  125. - Home page https://lemmasoftware.org
  126. - Git repository https://git.lemmasoftware.org
  127. - Broken Spoke Develpment http://numericalgeo.com
  128. - CGEM at the Coloroado School of Mines http://geophysics.mines.edu/cgem/
  129. - EGI at the Eniversity of Utah http://egi.utah.edu/
  130. **/
  131. #pragma once
  132. #ifndef __LEMMA_H
  133. #define __LEMMA_H
  134. #include <LemmaConfig.h>
  135. // Include some basic stuff that will always be needed
  136. #include <iostream>
  137. #include <iomanip>
  138. #include <complex>
  139. #include <fstream>
  140. #include <string>
  141. #include <vector>
  142. #include <stdexcept>
  143. #include <sstream>
  144. #include <Eigen/Core>
  145. #include <cstddef>
  146. #include <Eigen/StdVector>
  147. #include <Eigen/Sparse>
  148. #include <unsupported/Eigen/FFT>
  149. //#include <unsupported/Eigen/SparseExtra>
  150. #include <Eigen/Geometry>
  151. /** \brief The only namespace used by Lemma
  152. *
  153. * \details The rational behind this namespace is that built-in
  154. * types should be used wherever possible, but not
  155. * not built-in names. This allows for code that is better
  156. * enacsulated and easier to modify. The typedefs and constants
  157. * specified here are defined so that
  158. * precision/inplimentation can easily be changed.
  159. * All floating precision types should be typedefed in this file
  160. * and should not be used natively within any code.
  161. * Lemma uses
  162. * the Eigen Matrix/Vector/Linear Algebra library.
  163. * <http://eigen.tuxfamily.org> and a lot of the namespece typedefs
  164. * are specifying Eigen types.
  165. */
  166. namespace Lemma {
  167. /// Real defines precision for the whole API, default is double
  168. #ifdef LEMMA_SINGLE_PRECISION
  169. typedef float Real;
  170. #else // ----- LEMMA_SINGLE_PRECISION -----
  171. typedef double Real;
  172. #endif // ----- not LEMMA_SINGLE_PRECISION -----
  173. /// Complex version of Real.
  174. typedef std::complex<Real> Complex;
  175. /// A 3 component Eigen vector of Reals
  176. typedef Eigen::Matrix<Real, 3, 1> Vector3r;
  177. /// A 3 X Dynamic Component Eigen matrix of Reals
  178. typedef Eigen::Matrix<Real, 3, Eigen::Dynamic> Vector3Xr;
  179. /// Variable length Eigen vector of Reals
  180. typedef Eigen::Matrix<Real, Eigen::Dynamic, 1> VectorXr;
  181. /// Variable length Eigen vector of integers (int)
  182. typedef Eigen::Matrix<int, Eigen::Dynamic, 1> VectorXi;
  183. /// Variable length Eigen vector of Complexes
  184. typedef Eigen::Matrix<Complex, Eigen::Dynamic, 1> VectorXcr;
  185. /// A 3 Component Eigen vector of Complexes
  186. typedef Eigen::Matrix<Complex, 3, 1> Vector3cr;
  187. /// A 3 X Dynamic Component Eigen matrix of Complexes
  188. typedef Eigen::Matrix<Complex, 3, Eigen::Dynamic> Vector3Xcr;
  189. /// Variable length Eigen Matrix of Reals
  190. typedef Eigen::Matrix<Real, Eigen::Dynamic, Eigen::Dynamic> MatrixXr;
  191. /// Variable length Eigen Matrix of ints
  192. typedef Eigen::Matrix<int, Eigen::Dynamic, Eigen::Dynamic> MatrixXi;
  193. /// Variable length Eigen vector of Complexes
  194. typedef Eigen::Matrix<Complex, Eigen::Dynamic, Eigen::Dynamic> MatrixXcr;
  195. ////////////////////////////////////////
  196. // Constants used across the programmes
  197. /// Restating the obvious, this is pi
  198. const Real PI = 4.0*atan(1.0);
  199. /// Permitivity of Free Space
  200. //const Real EPSILON0 = 8.854187817e-12;
  201. const Real EPSILON0 = 8.854187817e-12;
  202. /// Permeability of free space
  203. const Real MU0 = 4.*PI*1e-7;
  204. /// 1/4 of \f$ \pi\f$
  205. const Real QPI = .25/PI;
  206. /// Some functions will convert units from SI (standard) to Gauss
  207. /// This is because NMR calculations are much more natural in Gauss
  208. enum MAGUNITS {TESLA, NANOTESLA, GAUSS};
  209. /// Unit of temperature entered
  210. enum TEMPUNITS {CELCIUS, KELVIN};
  211. /// Unit of time entered
  212. enum TIMEUNITS {SEC, MILLISEC, MICROSEC, NANOSEC, PICOSEC};
  213. /// Unit of time entered
  214. enum FREQUENCYUNITS {HZ, KHZ, MHZ, GHZ};
  215. /// FEM coil relative orientations
  216. enum FEMCOILORIENTATION {COAXIAL, COPLANAR};
  217. /// General orientation relative to coordinate system
  218. enum ORIENTATION {X, Y, Z, NX, NY, NZ};
  219. /// Type of field
  220. enum FIELDTYPE {HFIELDREAL, HFIELDIMAG, EFIELDREAL, EFIELDIMAG};
  221. /// Compenent of vector field
  222. enum FIELDCOMPONENT {XCOMPONENT=0, YCOMPONENT=1, ZCOMPONENT=2};
  223. /// Spatial component of vector
  224. enum SPATIALCOORDINANT {XCOORD=0, YCOORD=1, ZCOORD=2};
  225. /** Evaluation method for Hankel integrals.
  226. * ANDERSON801 Walt Anderson's 801 point filter
  227. * CHAVE Alan Chave's gaussian quadrature integration method
  228. * FHTKEY201 Key's 201 point filter
  229. * FHTKEY201 Key's 101 point filter
  230. * FHTKEY51 Key's 51 point filter
  231. * QWEKEY Key's Gaussian quadrature integration method
  232. */
  233. enum HANKELTRANSFORMTYPE { ANDERSON801, CHAVE, FHTKEY201, FHTKEY101, FHTKEY51, QWEKEY,
  234. FHTKONG61, FHTKONG121, FHTKONG241, IRONS };
  235. /** Enum is OK because these are the only physically possible sources.
  236. @param NOSOURCETYPE is default.
  237. @param ELECTRICDIPOLE is an electric dipole
  238. @param MAGNETICDIPOLE is a magnetic dipole
  239. */
  240. enum DipoleSourceType {NOSOURCETYPE, GROUNDEDELECTRICDIPOLE, UNGROUNDEDELECTRICDIPOLE, MAGNETICDIPOLE};
  241. /// Only three polarizations are supported. They may be summed to
  242. /// approximate others
  243. /// @param NOPOLARISATION is uninitialized, default value
  244. /// @param XPOLARISATION is a dipole oriented in the x direction
  245. /// @param YPOLARISATION is a dipole oriented in the y direction
  246. /// @param ZPOLARISATION is a dipole oriented in the z direction
  247. enum DipoleSourcePolarisation {NOPOLARISATION, XPOLARISATION,
  248. YPOLARISATION, ZPOLARISATION};
  249. /// The polarity may be either negative or positinve
  250. enum DipoleSourcePolarity {NEGATIVE, POSITIVE};
  251. /** The fields to make calculations on
  252. */
  253. enum FIELDCALCULATIONS {E, H, BOTH};
  254. /** Windowing function type
  255. */
  256. enum WINDOWTYPE { HAMMING, /*!< A hamming window */
  257. HANNING, /*!< A hanning window */
  258. RECTANGULAR /*!< Rectangular window */
  259. };
  260. }
  261. #endif // __Lemma_H