Page MenuHomec4science
Files F62067282

norm_utilities.py
No OneTemporary
Actions

File Metadata

Created
Fri, May 10, 17:21

norm_utilities.py
View Options

# Copyright (C) 2018 by the RROMPy authors
#
# This file is part of RROMPy.
#
# RROMPy is free software: you can redistribute it and/or modify
# it under the terms of the GNU Lesser General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# RROMPy is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Lesser General Public License for more details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with RROMPy. If not, see <http://www.gnu.org/licenses/>.
#
from abc import abstractmethod
import numpy as np
from copy import deepcopy as copy
from rrompy.utilities.base.types import Np1D, Np2D, DictAny
from rrompy.solver.linear_solver import setupSolver
from rrompy.utilities.exception_manager import RROMPyException
__all__ = ['Np2DLike', 'Np2DLikeEye', 'Np2DLikeInv', 'Np2DLikeInvLowRank',
'normEngine']
@abstractmethod
class Np2DLike:
def dot(self, u:Np2D) -> Np2D:
pass
class Np2DLikeEye(Np2DLike):
def dot(self, u:Np2D) -> Np2D:
return u
class Np2DLikeInv(Np2DLike):
def __init__(self, K:Np2D, M:Np2D, solverType:str, solverArgs:DictAny):
self.K, self.M, self.MH = K, M, M.T.conj()
self.solver, self.solverArgs = setupSolver(solverType, solverArgs)
def dot(self, u:Np2D) -> Np2D:
return self.MH.dot(self.solver(self.K, self.M.dot(u), self.solverArgs))
class Np2DLikeInvLowRank(Np2DLike):
def __init__(self, K:Np2D, M:Np2D, solverType:str, solverArgs:DictAny,
rank:int, oversampling : int = 10, seed : int = 420):
if rank > M.shape[1]:
raise RROMPyException(("Cannot select compressed rank larger than "
"original size."))
if oversampling < 0:
raise RROMPyException("Oversampling parameter must be positive.")
HF = Np2DLikeInv(K, M, solverType, solverArgs)
np.random.seed(seed)
xs = np.random.randn(M.shape[1], rank + oversampling)
samples = HF.dot(xs)
Q, _ = np.linalg.qr(samples, mode = "reduced")
R = HF.dot(Q).T.conj() # assuming HF (i.e. K) hermitian...
U, s, Vh = np.linalg.svd(R)
self.L = Q.dot(U[:, : rank]) * s[: rank]
self.R = Vh[: rank, :]
def dot(self, u:Np2D) -> Np2D:
return self.L.dot(self.R.dot(u))
class normEngine:
def __init__(self, energyNormMatrix:Np2D):
self.energyNormMatrix = copy(energyNormMatrix)
def innerProduct(self, u:Np2D, v:Np2D, onlyDiag : bool = False) -> Np2D:
if not isinstance(u, (np.ndarray,)): u = u.data
if not isinstance(v, (np.ndarray,)): v = v.data
if onlyDiag:
return np.sum(self.energyNormMatrix.dot(u) * v.conj(), axis = 0)
return v.T.conj().dot(self.energyNormMatrix.dot(u))
def norm(self, u:Np2D) -> Np1D:
return np.power(np.abs(self.innerProduct(u, u, onlyDiag = True)), .5)

Event Timeline

c4science · Help