# Source code for irfpy.util.graminan

''' Module for Gramian matrix <http://en.wikipedia.org/wiki/Gramian_matrix>_

'''
import logging
_logger = logging.getLogger(__name__)
import numpy as np

[docs]class Graminan:
def __init__(self, a_m):
r''' Instance granian.

:param a_m: M arrays with dimension of N.
Real space only considered.

Consider 3-D real space (N=3). Take 2 vector (M=2).
Then, Graminan can be obtained via

>>> a1 = np.array([1., 0, -2])
>>> a2 = np.array([0, -2., 1])
>>> G2 = Graminan([a1, a2])
>>> print(G2.m)
2
>>> print(G2.n)
3
>>> print(G2.matrix)
[[ 5. -2.]
[-2.  5.]]

Determinant of G highly relate to the area and volume.
For this case, only two vector povides only area.

.. math::

S = \sqrt{det(G_2)}

For m=3 case, you can get volume as

.. math::

V = \sqrt{det(G_3)}

Probably for m=4 case, you can get super volume as

.. math::

V_4 = \sqrt{det(G_4)}

These functionalities are implemented as :meth:get_volume function.
'''

self.m = len(a_m)
self.n = len(a_m[0])

if self.m > self.n:
_logger.error('The Graminan should have m<=n, while m=%d n=%d' % (self.m, self.n))
raise ValueError('m=%d n=%d, should be m<=n' % (self.m, self.n))

m = self.m

self.matrix = np.zeros([m, m])
for ix in range(m):
for iy in range(m):
self.matrix[ix, iy] = np.dot(a_m[ix], a_m[iy])

[docs]    def get_volume(self):
''' Return the volume (generalized).

For m=2 case, area (of trapezoid) is returned.
For m=3 case, volume (hexahedron) is returend.
Probably for m=4 case, super volume is returend.
'''
return np.sqrt(np.linalg.det(self.matrix))

import unittest
import doctest

[docs]def doctests():
return unittest.TestSuite((
doctest.DocTestSuite(),
))
if __name__ == '__main__':
unittest.main(defaultTest='doctests')