Source code for seapy.system

"""
System
======

The System class is the main class and contains methods for solving the SEA model.

.. autoclass:: seapy.system.System

.. autoclass:: seapy.system.Frequency
.. autoclass:: seapy.system.Band

"""

import math
import cmath
import numpy as np

import warnings
import logging
from weakref import WeakSet
import weakref

from .base import Spectrum
from seapy.objects_map import objects_map

from tabulate import tabulate
from toolz import count
import pandas as pd
from .tools import plot


[docs]class System(object): """ The System class contains methods for solving the model. """ SORT = 'System' reference_pressure = 2.0e-5 """Reference pressure :math:`p_{0}`. """ reference_velocity = 5.0e-8 """Reference velocity :math:`v_{0}`. """ reference_power = 1.0e-12 """Reference power :math:`P_{0}`. """ solved = False """ Switch indicating whether the system (modal energies) were solved or not. """ @property
[docs] def objects(self): """All objects in the system. :returns: Generator of objects. :rtype: :class:`types.GeneratorType` """ yield from (self.getObject(obj.name) for obj in self._objects)
@property
[docs] def components(self): """All components in the system. :returns: Generator of components. :rtype: :class:`python.types.GeneratorType` """ yield from (obj for obj in self.objects if obj.SORT=='Component')
@property
[docs] def subsystems(self): """All subsystems in the system. :returns: Generator of subsystems. :rtype: :class:`python.types.GeneratorType` """ yield from (obj for obj in self.objects if obj.SORT=='Subsystem')
@property
[docs] def junctions(self): """All junctions in the system. :returns: Generator of junctions. :rtype: :class:`python.types.GeneratorType` """ yield from (obj for obj in self.objects if obj.SORT=='Junction')
@property
[docs] def couplings(self): """All couplings in the system. :returns: Generator of couplings. :rtype: :class:`python.types.GeneratorType` """ yield from (obj for obj in self.objects if obj.SORT=='Coupling')
@property
[docs] def materials(self): """All materials in the system. :returns: Generator of materials. :rtype: :class:`python.types.GeneratorType` """ yield from (obj for obj in self.objects if obj.SORT=='Material')
@property
[docs] def excitations(self): """All excitations in the system. :returns: Generator of excitations. :rtype: :class:`python.types.GeneratorType` """ yield from (obj for obj in self.objects if obj.SORT=='Excitation')
[docs] def __init__(self): """Constructor. """ self.frequency = Frequency(weakref.proxy(self)) """ Frequency object """ self._objects = list() """Private set of objects this SEA model consists of. """
def __del__(self): for obj in self.objects: self.removeObject(obj.name) def _getRealObject(self, name): """ Get real object by name. :param name: Name of `object`. :returns: Real `object`. """ name = name if isinstance(name, str) else name.name for obj in self._objects: if name == obj.name: return obj else: raise ValueError("Unknown name. Cannot get object.")
[docs] def getObject(self, name): """ Get object by name. :param name: Name of `object`. :returns: Proxy to `object`. """ name = name if isinstance(name, str) else name.name for obj in self._objects: if name == obj.name: return weakref.proxy(obj) else: raise ValueError("Unknown name. Cannot get object.")
[docs] def removeObject(self, name): """ Delete object from SEA model. :param name: Name of `object`. :returns: Proxy to `object`. """ obj = obj._getRealObject(name) for obj in self._objects: if name == obj.name: self._objects.remove(obj)
def _addObject(self, name, model, **properties): """Add object to SEA model. :param name: Name of `object`. :param model: Model or type of `object`. :param properties: Other properties specific to `object`. """ #if name not in self._objects: try: obj = model(name, weakref.proxy(self), **properties) # Add hidden hard reference except KeyError: raise ValueError("Model does not exist. Cannot create object.") self._objects.append(obj) return self.getObject(obj.name)
[docs] def addComponent(self, name, model, **properties): """Add component to SEA model. :param name: Name of component. :param model: Model or type of component. See :attr:`seapy.components.components_map`. :param properties: Other properties specific to the component. """ obj = self._addObject(name, objects_map['component'][model] , **properties) obj._addSubsystems() return obj
[docs] def addJunction(self, name, model, **properties): """Add junction to SEA model.""" obj = self._addObject(name, objects_map['junction'][model], **properties) #obj._updateCouplings() return obj
[docs] def addMaterial(self, name, model, **properties): """Add material to SEA model.""" obj = self._addObject(name, objects_map['material'][model], **properties) return obj
[docs] def addCoupling(self, name, model, **properties): """Add material to SEA model.""" obj = self._addObject(name, objects_map['coupling'][model], **properties) return obj
[docs] def power_balance_matrix(self): """Power balance matrix as function of frequency. :param subsystems: is a list of subsystems. Reason to give the list as argument instead of using self.subsystems is that that list might change during execution. :type subsystems: list :param f: is the index of the center frequency of the frequency band :type f: int :rtype: :class:`numpy.ndarray` See Craik, equation 6.21, page 155. Yields power balance matrices. """ #logging.info("Creating matrix for centerfrequency {}".format(self.frequency.center[f])) subsystems = [subsystem for subsystem in self.subsystems if subsystem.included==True] for f in range(self.frequency.amount): B = np.zeros((len(subsystems), len(subsystems)), dtype=float) for j, subsystem_j in enumerate(subsystems): # Row j for i, subsystem_i in enumerate(subsystems): # Column i loss_factor = 0.0 if i==j: loss_factor = + subsystem_i.tlf[f] # Total loss factor. else: ####Take the coupling loss factor from subsystem i to subsystem j. Negative x = list(set(subsystem_i.linked_couplings_from).intersection(set(subsystem_j.linked_couplings_to))) ##if not x: ## Use the relation consistency relationship? ##pass ##print 'error. No coupling?' #print(len(x)) if len(x)==1: #print("We have one!") #print(len(x)) loss_factor = - x[0].clf[f] #coupling = x[0] ##loss_factor = - coupling.clf[f] del x B[j,i] = loss_factor * subsystem_i.modal_density[f] i+=1 j+=1 logging.info('Matrix created.') logging.info(B) yield B
[docs] def power_vector(self): """Vector of input power normalized with angular frequency. See Craik, equation 6.21, page 155 Yields power input vectors. """ subsystems = [subsystem for subsystem in self.subsystems if subsystem.included==True] for f in range(self.frequency.amount): power_input = [subsystem.power_input[f] / self.frequency.angular[f] for subsystem in subsystems] yield np.array(power_input)
[docs] def clearResults(self): """Clear the results. Reset modal energies. Set :attr:`solved` to False. :rtype: None """ logging.info('Clearing results...') for subsystem in self.subsystems: del subsystem.modal_energy self.solved = False logging.info('Cleared results.')
[docs] def solveSystem(self): # Put the actual solving in a separate thread? """Solve modal powers. :rtype: :func:`bool` This method solves the modal energies for every subsystem. The method :meth:`createMatrix` is called for every frequency band to construct a matrix of :term:`loss factors` and :term:`modal densities`. """ logging.info('Solving system...') self.clean() subsystems = [subsystem for subsystem in self.subsystems if subsystem.included==True] for f, (p, B) in enumerate(zip(self.power_vector(), self.power_balance_matrix())): if self.frequency.enabled[f]: #try: modal_energy = np.linalg.solve(B, p) # Left division results in the modal energies. #except np.linalg.linalg.LinAlgError as e: # If there is an error solving the matrix, then quit right away. #warnings.warn( repr(e) ) #return False for i, subsystem in enumerate(subsystems): subsystem.modal_energy[f] = modal_energy[i] #del modal_energy, power_input, LF self.solved = True logging.info('System solved.') return True
[docs] def clean(self): """Reset modal energy to zero in all subsystems. """ for subsystem in self.subsystems: subsystem.modal_energy = np.zeros(self.frequency.amount) #def info(self, sort=None, fields=None): #"""Print information about objects of type sort. By default all types are returned."""
def plot(self, objects, quantity, yscale='linear'): items = (getattr(obj, quantity) for obj in (self.getObject(o) for o in objects)) return plot(self.frequency.center, items, quantity, self.frequency.enabled, yscale=yscale)
[docs] def info(self, objects, attribute, tablefmt='simple'): """Print single attribute of objects. """ objects = (self.getObject(obj) for obj in objects) data = {obj.name: getattr(obj, attribute) for obj in objects if hasattr(obj, attribute)} df = pd.DataFrame(data, index=self.frequency.center.astype('int')).T return df
[docs] def info_objects(self): """Print information about objects. """ header = ['Name', 'Included', 'Enabled', 'Class', 'Sort'] data = ((obj.name, obj.included, obj.enabled, obj.__class__.__name__, obj.SORT) for obj in self.objects) return tabulate(data, headers=header, tablefmt=tablefmt)
[docs] def info_materials(self): """Print information about materials. """ header = ['Name', 'Included', 'Enabled', 'Class', 'Density', 'Components'] data = ((obj.name, obj.included, obj.enabled, obj.__class__.__name__, obj.density, count(obj.linked_components)) for obj in self.materials) return tabulate(data, headers=header, tablefmt=tablefmt)
[docs] def components_info(self, tablefmt='simple'): """Print information about components. """ header = ['Name', 'Included', 'Enabled', 'Class', 'Material', 'Volume', 'Mass', 'Subsystems', 'Junctions'] data = ((obj.name, obj.included, obj.enabled, obj.__class__.__name__, obj.material.name, obj.volume, obj.mass, count(obj.linked_subsystems), count(obj.linked_junctions)) for obj in self.components) return tabulate(data, headers=header, tablefmt=tablefmt)
[docs] def subsystems_info(self, tablefmt='simple'): """Print information about subsystems. """ header = ['Name', 'Included', 'Enabled', 'Class', 'Component', 'Couplings - From', 'Couplings - To', 'Excitations'] data = ((obj.name, obj.included, obj.enabled, obj.__class__.__name__, obj.component.name, count(obj.linked_couplings_from), count(obj.linked_couplings_to), count(obj.linked_excitations)) for obj in self.subsystems) return tabulate(data, headers=header, tablefmt=tablefmt)
[docs] def excitations_info(self, tablefmt='simple'): """Print information about excitations. """ header = ['Name', 'Included', 'Enabled', 'Class', 'Subsystem'] data = ((obj.name, obj.included, obj.enabled, obj.__class__.__name__, obj.subsystem.name) for obj in self.excitations) return tabulate(data, headers=header, tablefmt=tablefmt)
[docs] def junctions_info(self, tablefmt='simple'): """Print information about junctions. """ header = ['Name', 'Included', 'Enabled', 'Class', 'Shape', 'Components', 'Couplings', 'Subsystems'] data = ((obj.name, obj.included, obj.enabled, obj.__class__.__name__, obj.shape, count(obj.components), count(obj.linked_couplings), count(obj.subsystems)) for obj in self.junctions) return tabulate(data, headers=header, tablefmt=tablefmt) #def couplings_info(self, tablefmt='simple'): #"""Print information about junctions. #""" #header = ['Name', 'included', 'Enabled', 'Class', 'Shape', 'Components', 'Couplings', 'Subsystems'] #data = ((obj.name, obj.included, obj.impedance, obj.__class__.__name__, obj.shape, count(obj.components), count(obj.linked_couplings), count(obj.subsystems)) for obj in self.excitations) #return tabulate(data, headers=header, tablefmt=tablefmt)
[docs]class Band(object): """Frequency band class."""
[docs] def __init__(self, lower=0.0, center=0.0, upper=0.0, enabled=False): self.lower = lower self.center = center self.upper = upper self.enabled = enabled
@property def bandwidth(self): return self.upper - self.lower @property def angular(self): return 2.0 * np.pi * self.center #class SpectrumDescriptorFrequency(object): #""" #""" #def __get__(self): #pass #def __set__(self): #pass
class Frequency(object): """New-style spectrum class.""" def __init__(self, system): self.system = system self._bands = list() def _spectrum(name): """Property to access the frequency bands as/using arrays.""" @property def prop(self): return np.array([getattr(band, name) for band in self._bands]) @prop.setter def prop(self, x): if len(x) == len(self._bands): """ When the given array has the same amount of items as there are frequency bands, we will fit them one on one. """ for new, band in zip(x, self._bands): setattr(band, name, new) else: """ If not, we will delete the old frequency bands, and create new ones.""" self._bands = list() for i in x: band = Band() setattr(band, name, i) self._bands.append(band) # Use self.addBand() instead!! return prop lower = _spectrum('lower') center = _spectrum('center') upper = _spectrum('upper') enabled = _spectrum('enabled') @property def bandwidth(self): return np.array([band.bandwidth for band in self._bands]) @property def angular(self): return np.array([band.angular for band in self._bands]) @property def amount(self): return len(self._bands) @property def spectra(self): """Generator to obtain all spectra in use in the SEA model.""" for obj in self.system.objects: for cls in obj.__class__.__mro__: for key, value in cls.__dict__.items(): if isinstance(value, Spectrum): yield (obj, key) def appendBand(self, **kwargs): """ Append frequency band. """ self.addBand(len(self._bands), **kwargs) def addBand(self, pos, **kwargs): """ Add frequency band. Inform all spectra of the change!! """ self._bands.insert(pos, Band(**kwargs)) # Create a new frequency band default = 0.0 # default value of array cells for obj, attr in self.iterspectra(): # Add a band to all spectra setattr(obj, attr, np.insert(getattr(obj, attr), pos, default)) def removeBand(self, pos): """ Remove frequency band. Inform all spectra of the change!! """ self._bands.pop(pos) # Remove the frequency band for obj, attr in self.iterspectra(): # Remove a band from all spectra setattr(obj, attr, np.delete(getattr(obj, attr), pos)) #@property #def upper(self): #return np.array([band.upper for band in self._bands]) #@property #def lower(self): #return np.array([band.lower for band in self._bands]) #@property #def center(self): #return np.array([band.center for band in self._bands]) #@center.setter #def center(self, x): #if len(x) == len(self.bands): #for new, band in zip(x, self.bands): #band.center = new ##else: #"""Add the required amount of bands.""" #@property #def enabled(self): #return np.array([band.enabled for band in self._bands]) #class Frequency(object): #""" #Abstract base class for handling different frequency settings. #""" #def __init__(self, system): #self._system = system #def _set_band(self, x, sort): #for i in ['_lower', '_center', '_upper']: #if len(getattr(self, i)) != len(x): #setattr(self, i, np.zeros(len(x))) #if len(getattr(self, '_enabled')) != len(x): #setattr(self, '_enabled', np.zeros(len(x), dtype=bool)) #setattr(self, sort, np.array(x)) #def _get_center(self): #return self._center #def _set_center(self, x): #self._set_band(x, '_center') #_center = np.array([0.0]) #center = property(fget=_get_center, fset=_set_center) #""" #Center frequencies of frequency bands. #""" #def _get_upper(self): #return self._upper #def _set_upper(self, x): #self._set_band(x, '_upper') #_upper = np.array([0.0]) #upper = property(fget=_get_upper, fset=_set_upper) #""" #Upper limit frequencies of frequency bands. #""" #def _set_lower(self, x): #self._set_band(x, '_lower') #def _get_lower(self): #return self._lower #_lower = np.array([0.0]) #lower = property(fget=_get_lower, fset=_set_lower) #""" #Lower limit frequencies of frequency bands. #""" #def _set_enabled(self, x): #self._set_band(x, '_enabled') #def _get_enabled(self): #return self._enabled #_enabled = np.array([False]) #enabled = property(fget=_get_upper, fset=_set_upper) #""" #Enabled frequency bands. #Modal powers will not be solved for disabled frequency bands. #""" #@property #def bandwidth(self): #"""Bandwidth of frequency bands, #:rtype: :class:`numpy.ndarray` #""" #return self.upper - self.lower #@property #def angular(self): #"""Angular frequency #:rtype: :class:`numpy.ndarray` #""" #return self.center * 2.0 * np.pi #@property #def amount(self): #"""Amount of frequency bands #:rtype: :func:`int` #""" #try: #return len(self.center) #except TypeError: #return 0