Source code for openmdao.lib.datatypes.domain.metrics

"""
Predefined classes used by :meth:`mesh_probe`.
Some typical fluid flow metric calculation classes are included here.
Others may be added via :meth:`register_metric`. A list of available metrics
is returned by :meth:`list_metrics`. Note that the name a metric is registered
under need not match its class name.

Metrics may be used with 1D, 2D, or 3D Cartesian coordinates. They may also
be used with polar (2D) or cylindrical (3D) coordinates. :meth:`calculate`
should be prepared for this.
"""

from math import sqrt

from openmdao.units.units import PhysicalQuantity

from openmdao.lib.datatypes.domain.zone import CYLINDRICAL

# Dictionary for metric information.
# Populated as metric classes are registered.
_METRICS = {}

# Legal `geometry` values.
_GEOMETRIES = ('volume', 'surface', 'curve', 'any')


[docs]def register_metric(name, cls, integrate, geometry='any'): """ Register a metric calculation class. name: string The name of the metric calculated. cls: class Class used to calculate metric. It must contain :meth:`calculate` and :meth:`dimensionalize`. Constructor arguments are `(zone, zone_name, reference_state)`. :meth:`calculate` will be called with `(loc, geom)`. `loc` contains indices into the zone variable arrays and `geom` is either the cell volume, a non-dimensional vector normal to the cell face with magnitude equal to its area, or the edge length; depending upon the type of region (volume, surface, or curve). :meth:`dimensionalize` is called with the accumulated value. It should return a :class:`PhysicalQuantity` for the dimensionalized value. integrate: bool If True, then calculated values are integrated, not averaged. geometry: string ('volume', 'surface', 'curve', or 'any') Specifies what type of geometry this metric is restricted to, if any. """ if geometry not in _GEOMETRIES: raise ValueError('geom must be one of %s' % _GEOMETRIES) _METRICS[name] = (cls, integrate, geometry)
[docs]def get_metric(name): """ Return ``(cls, integrate, geometry)``. """ return _METRICS[name]
[docs]def list_metrics(): """ Return list of valid metric names. """ return sorted(_METRICS.keys())
[docs]def create_scalar_metric(var_name): """ Creates a minimal metric calculation class for `var_name` and registers it. This can be used for averaging scalar flow variables without having to statically create a class and register it. The created class does not support dimensionalization. """ cls_name = var_name.capitalize() exec ''' class %(cls_name)s(object): """ Computes %(var_name)s. """ def __init__(self, zone, zone_name, reference_state): self.%(var_name)s = zone.flow_solution.%(var_name)s.item def calculate(self, loc, length): """ Return metric value. """ return self.%(var_name)s(*loc) def dimensionalize(self, value): """ Return dimensional `value`. """ raise NotImplementedError('Dimensional %(var_name)s') register_metric('%(var_name)s', %(cls_name)s, False, 'any') ''' % {'var_name': var_name, 'cls_name': cls_name}
[docs]class Area(object): """ Computes area of mesh surface. """ def __init__(self, zone, zone_name, reference_state): if reference_state is None: self.aref = 1. else: try: lref = reference_state['length_reference'] except KeyError: raise AttributeError("For area, reference_state is missing" " 'length_reference'.") aref = lref * lref self.units = aref.get_unit_name() self.aref = aref.value
[docs] def calculate(self, loc, normal): """ Return metric value. """ sc1, sc2, sc3 = normal sc1 *= self.aref sc2 *= self.aref sc3 *= self.aref return sqrt(sc1*sc1 + sc2*sc2 + sc3*sc3)
[docs] def dimensionalize(self, value): """ Return dimensional `value`. """ return PhysicalQuantity(value, self.units)
register_metric('area', Area, True, 'surface')
[docs]class Length(object): """ Computes length of mesh curve. """ def __init__(self, zone, zone_name, reference_state): if reference_state is None: self.units = None else: try: lref = reference_state['length_reference'] except KeyError: raise AttributeError("For length, reference_state is missing" " 'length_reference'.") self.units = lref.get_unit_name() self.lref = lref.value
[docs] def calculate(self, loc, length): """ Return metric value. """ return length * self.lref
[docs] def dimensionalize(self, value): """ Return dimensional `value`. """ return PhysicalQuantity(value, self.units)
register_metric('length', Length, True, 'curve')
[docs]class MassFlow(object): """ Computes mass flow across a mesh surface. """ def __init__(self, zone, zone_name, reference_state): flow = zone.flow_solution cylindrical = zone.coordinate_system == CYLINDRICAL try: momentum = flow.momentum except AttributeError: raise AttributeError("For mass_flow, zone %s is missing 'momentum'." % zone_name) if reference_state is None: self.aref = 1. self.momref = 1. self.wref = None else: try: lref = reference_state['length_reference'] pref = reference_state['pressure_reference'] rgas = reference_state['ideal_gas_constant'] tref = reference_state['temperature_reference'] except KeyError: vals = ('length_reference', 'pressure_reference', 'ideal_gas_constant', 'temperature_reference') raise AttributeError('For mass_flow, reference_state is missing' ' one or more of %s.' % (vals,)) rhoref = pref / rgas / tref vref = (rgas * tref).sqrt() momref = rhoref * vref self.wref = momref * lref * lref aref = lref * lref self.aref = aref.value self.momref = momref.value if cylindrical: self.mom_c1 = None if momentum.z is None else momentum.z.item self.mom_c2 = momentum.r.item self.mom_c3 = momentum.t.item else: self.mom_c1 = momentum.x.item self.mom_c2 = None if momentum.y is None else momentum.y.item self.mom_c3 = None if momentum.z is None else momentum.z.item
[docs] def calculate(self, loc, normal): """ Return metric value. """ rvu = 0. if self.mom_c1 is None else self.mom_c1(*loc) * self.momref rvv = 0. if self.mom_c2 is None else self.mom_c2(*loc) * self.momref rvw = 0. if self.mom_c3 is None else self.mom_c3(*loc) * self.momref sc1, sc2, sc3 = normal sc1 *= self.aref sc2 *= self.aref sc3 *= self.aref return rvu*sc1 + rvv*sc2 + rvw*sc3
[docs] def dimensionalize(self, value): """ Dimensionalize `value`. """ return PhysicalQuantity(value, self.wref.get_unit_name())
register_metric('mass_flow', MassFlow, True, 'surface')
[docs]class CorrectedMassFlow(object): """ Computes corrected mass flow across a mesh surface. """ def __init__(self, zone, zone_name, reference_state): flow = zone.flow_solution cylindrical = zone.coordinate_system == CYLINDRICAL # 'pressure' required until we can determine dimensionalized # static pressure from 'Q' variables. try: self.density = flow.density.item momentum = flow.momentum self.pressure = flow.pressure.item except AttributeError: vnames = ('density', 'momentum', 'pressure') raise AttributeError('For corrected_mass_flow, zone %s is missing' ' one or more of %s.' % (zone_name, vnames)) try: self.gam = flow.gamma.item except AttributeError: self.gam = None # Use passed-in scalar gamma. if reference_state is None: raise ValueError('corrected_mass_flow must have units specified') try: lref = reference_state['length_reference'] pref = reference_state['pressure_reference'] rgas = reference_state['ideal_gas_constant'] tref = reference_state['temperature_reference'] if self.gam is None: self.gamma = reference_state['specific_heat_ratio'].value except KeyError: vals = ('length_reference', 'pressure_reference', 'ideal_gas_constant', 'temperature_reference', 'specific_heat_ratio') raise AttributeError('For corrected_mass_flow, reference_state is' ' missing one or more of %s.' % (vals,)) rhoref = pref / rgas / tref vref = (rgas * tref).sqrt() momref = rhoref * vref self.wref = momref * lref * lref pstd = PhysicalQuantity(14.696, 'psi') pstd.convert_to_unit(pref.get_unit_name()) tstd = PhysicalQuantity(518.67, 'degR') tstd.convert_to_unit(tref.get_unit_name()) aref = lref * lref self.aref = aref.value self.pref = pref.value self.rgas = rgas.value self.rhoref = rhoref.value self.momref = momref.value self.pstd = pstd.value self.tstd = tstd.value if cylindrical: self.mom_c1 = None if momentum.z is None else momentum.z.item self.mom_c2 = momentum.r.item self.mom_c3 = momentum.t.item else: self.mom_c1 = momentum.x.item self.mom_c2 = None if momentum.y is None else momentum.y.item self.mom_c3 = None if momentum.z is None else momentum.z.item
[docs] def calculate(self, loc, normal): """ Return metric value. """ rho = self.density(*loc) * self.rhoref rvu = 0. if self.mom_c1 is None else self.mom_c1(*loc) * self.momref rvv = 0. if self.mom_c2 is None else self.mom_c2(*loc) * self.momref rvw = 0. if self.mom_c3 is None else self.mom_c3(*loc) * self.momref ps = self.pressure(*loc) * self.pref if self.gam is not None: gamma = self.gam(*loc) else: gamma = self.gamma sc1, sc2, sc3 = normal sc1 *= self.aref sc2 *= self.aref sc3 *= self.aref w = rvu*sc1 + rvv*sc2 + rvw*sc3 u2 = (rvu*rvu + rvv*rvv + rvw*rvw) / (rho*rho) a2 = (gamma * ps) / rho mach2 = u2 / a2 ts = ps / (rho * self.rgas) tt = ts * (1. + (gamma-1.)/2. * mach2) pt = ps * pow(1. + (gamma-1.)/2. * mach2, gamma/(gamma-1.)) return w * sqrt(tt/self.tstd) / (pt/self.pstd)
[docs] def dimensionalize(self, value): """ Dimensionalize `value`. """ return PhysicalQuantity(value, self.wref.get_unit_name())
register_metric('corrected_mass_flow', CorrectedMassFlow, True, 'surface')
[docs]class StaticPressure(object): """ Computes weighted static pressure for a mesh region. """ def __init__(self, zone, zone_name, reference_state): flow = zone.flow_solution cylindrical = zone.coordinate_system == CYLINDRICAL try: # Some codes have this directly available. self.pressure = flow.pressure.item except AttributeError: self.pressure = None try: # Look for typical Q variables. self.density = flow.density.item momentum = flow.momentum self.energy = flow.energy_stagnation_density.item except AttributeError: vnames = ('pressure', 'density', 'momentum', 'energy_stagnation_density') raise AttributeError('For pressure, zone %s is missing' ' one or more of %s.' % (zone_name, vnames)) try: self.gam = flow.gamma.item except AttributeError: self.gam = None # Use passed-in scalar gamma. if reference_state is None: self.pref = 1. self.rhoref = 1. self.momref = 1. self.e0ref = 1. self.units = None else: if self.pressure is not None: try: pref = reference_state['pressure_reference'] except KeyError: raise AttributeError("For pressure, reference_state is" " missing 'pressure_reference'.") else: raise NotImplementedError('Get dimensional pressure from' ' Q variables') # pref = reference_state['pressure_reference'] # rgas = reference_state['ideal_gas_constant'] # tref = reference_state['temperature_reference'] # if self.gam is None: # self.gamma = reference_state['specific_heat_ratio'].value # except AttributeError: # vnames = ('ideal_gas_constant', 'pressure_reference', # 'temperature_reference', 'specific_heat_ratio') # raise AttributeError('For static pressure, reference_state is' # ' missing one or more of %s.' % vnames) # rhoref = pref / rgas / tref # vref = (rgas * tref).sqrt() # momref = rhoref * vref # self.rhoref = rhoref.value # self.momref = momref.value # self.e0ref = ? self.units = pref.get_unit_name() self.pref = pref.value if self.pressure is None: if cylindrical: self.mom_c1 = None if momentum.z is None else momentum.z.item self.mom_c2 = momentum.r.item self.mom_c3 = momentum.t.item else: self.mom_c1 = momentum.x.item self.mom_c2 = None if momentum.y is None else momentum.y.item self.mom_c3 = None if momentum.z is None else momentum.z.item
[docs] def calculate(self, loc, geom): """ Return metric value. """ if self.pressure is not None: return self.pressure(*loc) * self.pref else: rho = self.density(*loc) * self.rhoref vu = 0. if self.mom_c1 is None else self.mom_c1(*loc) * self.momref / rho vv = 0. if self.mom_c2 is None else self.mom_c2(*loc) * self.momref / rho vw = 0. if self.mom_c3 is None else self.mom_c3(*loc) * self.momref / rho e0 = self.energy(*loc) * self.e0ref / rho if self.gam is not None: gamma = self.gam(*loc) else: gamma = self.gamma return (gamma-1.) * rho * (e0 - 0.5*(vu*vu + vv*vv + vw*vw))
[docs] def dimensionalize(self, value): """ Dimensionalize `value`. """ return PhysicalQuantity(value, self.units)
register_metric('pressure', StaticPressure, False)
[docs]class TotalPressure(object): """ Computes weighted total pressure for a mesh region. """ def __init__(self, zone, zone_name, reference_state): flow = zone.flow_solution cylindrical = zone.coordinate_system == CYLINDRICAL try: self.density = flow.density.item momentum = flow.momentum except AttributeError: vnames = ('density', 'momentum') raise AttributeError('For pressure_stagnation, zone %s is missing' ' one or more of %s.' % (zone_name, vnames)) try: self.pressure = flow.pressure.item except AttributeError: self.pressure = None try: self.energy = flow.energy_stagnation_density.item except AttributeError: vnames = ('pressure', 'energy_stagnation_density') raise AttributeError('For pressure_stagnation, zone %s is missing' ' one or more of %s.' % (zone_name, vnames)) try: self.gam = flow.gamma.item except AttributeError: self.gam = None # Use passed-in scalar gamma. if reference_state is None: self.pref = 1. self.rhoref = 1. self.momref = 1. self.e0ref = 1. self.units = None else: try: pref = reference_state['pressure_reference'] rgas = reference_state['ideal_gas_constant'] tref = reference_state['temperature_reference'] if self.gam is None: self.gamma = reference_state['specific_heat_ratio'].value except KeyError: vals = ('pressure_reference', 'ideal_gas_constant', 'temperature_reference', 'specific_heat_ratio') raise AttributeError('For pressure_stagnation, reference_state' ' is missing one or more of %s.' % (vals,)) if self.pressure is None: raise NotImplementedError('Get dimensional pressure_stagnation from' ' Q variables') rhoref = pref / rgas / tref vref = (rgas * tref).sqrt() momref = rhoref * vref self.rhoref = rhoref.value self.momref = momref.value # self.e0ref = ? self.units = pref.get_unit_name() self.pref = pref.value if cylindrical: self.mom_c1 = None if momentum.z is None else momentum.z.item self.mom_c2 = momentum.r.item self.mom_c3 = momentum.t.item else: self.mom_c1 = momentum.x.item self.mom_c2 = None if momentum.y is None else momentum.y.item self.mom_c3 = None if momentum.z is None else momentum.z.item
[docs] def calculate(self, loc, geom): """ Return metric value. """ rho = self.density(*loc) * self.rhoref vu = 0. if self.mom_c1 is None else self.mom_c1(*loc) * self.momref / rho vv = 0. if self.mom_c2 is None else self.mom_c2(*loc) * self.momref / rho vw = 0. if self.mom_c3 is None else self.mom_c3(*loc) * self.momref / rho if self.gam is not None: gamma = self.gam(*loc) else: gamma = self.gamma u2 = vu*vu + vv*vv + vw*vw if self.pressure is not None: ps = self.pressure(*loc) * self.pref else: e0 = self.energy(*loc) * self.e0ref / rho ps = (gamma-1.) * rho * (e0 - 0.5*u2) a2 = (gamma * ps) / rho mach2 = u2 / a2 return ps * pow(1. + (gamma-1.)/2. * mach2, gamma/(gamma-1.))
[docs] def dimensionalize(self, value): """ Dimensionalize `value`. """ return PhysicalQuantity(value, self.units)
register_metric('pressure_stagnation', TotalPressure, False)
[docs]class StaticTemperature(object): """ Computes weighted static temperature for a mesh region. """ def __init__(self, zone, zone_name, reference_state): flow = zone.flow_solution cylindrical = zone.coordinate_system == CYLINDRICAL try: self.density = flow.density.item except AttributeError: raise AttributeError('For temperature, zone %s is missing' ' density.' % zone_name) try: self.pressure = flow.pressure.item except AttributeError: self.pressure = None try: # Look for typical Q variables. momentum = flow.momentum self.energy = flow.energy_stagnation_density.item except AttributeError: vnames = ('pressure', 'momentum', 'energy_stagnation_density') raise AttributeError('For temperature, zone %s is missing' ' one or more of %s.' % (zone_name, vnames)) try: self.gam = flow.gamma.item except AttributeError: self.gam = None # Use passed-in scalar gamma. if reference_state is None: self.pref = 1. self.rhoref = 1. self.momref = 1. self.e0ref = 1. self.rgas = 1. self.gamma = 1.4 self.tref = None else: try: pref = reference_state['pressure_reference'] rgas = reference_state['ideal_gas_constant'] tref = reference_state['temperature_reference'] except KeyError: vals = ('pressure_reference', 'ideal_gas_constant', 'temperature_reference') raise AttributeError('For temperature, reference_state is missing' ' one or more of %s.' % (vals,)) if self.pressure is None: if self.gam is None: try: self.gamma = reference_state['specific_heat_ratio'].value except KeyError: raise AttributeError('For temperature, reference_state is' ' missing specific_heat_ratio') raise NotImplementedError('Get dimensional temperature from' ' Q variables') rhoref = pref / rgas / tref self.pref = pref.value self.rgas = rgas.value self.rhoref = rhoref.value # self.e0ref = ? self.tref = tref if self.pressure is None: if cylindrical: self.mom_c1 = None if momentum.z is None else momentum.z.item self.mom_c2 = momentum.r.item self.mom_c3 = momentum.t.item else: self.mom_c1 = momentum.x.item self.mom_c2 = None if momentum.y is None else momentum.y.item self.mom_c3 = None if momentum.z is None else momentum.z.item
[docs] def calculate(self, loc, geom): """ Return metric value. """ rho = self.density(*loc) * self.rhoref if self.pressure is not None: ps = self.pressure(*loc) * self.pref else: vu = 0. if self.mom_c1 is None else self.mom_c1(*loc) * self.momref / rho vv = 0. if self.mom_c2 is None else self.mom_c2(*loc) * self.momref / rho vw = 0. if self.mom_c3 is None else self.mom_c3(*loc) * self.momref / rho e0 = self.energy(*loc) * self.e0ref / rho if self.gam is not None: gamma = self.gam(*loc) else: gamma = self.gamma ps = (gamma-1.) * rho * (e0 - 0.5*(vu*vu + vv*vv + vw*vw)) return ps / (rho * self.rgas)
[docs] def dimensionalize(self, value): """ Dimensionalize `value`. """ return PhysicalQuantity(value, self.tref.get_unit_name())
register_metric('temperature', StaticTemperature, False)
[docs]class TotalTemperature(object): """ Computes weighted total temperature for a mesh region. """ def __init__(self, zone, zone_name, reference_state): flow = zone.flow_solution cylindrical = zone.coordinate_system == CYLINDRICAL try: self.density = flow.density.item momentum = flow.momentum except AttributeError: vnames = ('density', 'momentum') raise AttributeError('For temperature_stagnation, zone %s is missing' ' one or more of %s.' % (zone_name, vnames)) try: self.pressure = flow.pressure.item except AttributeError: self.pressure = None try: self.energy = flow.energy_stagnation_density.item except AttributeError: vnames = ('pressure', 'energy_stagnation_density') raise AttributeError('For temperature_stagnation, zone %s is' ' one or more of %s.' % (zone_name, vnames)) try: self.gam = flow.gamma.item except AttributeError: self.gam = None # Use passed-in scalar gamma. if reference_state is None: self.pref = 1. self.rgas = 1. self.rhoref = 1. self.momref = 1. self.e0ref = 1. else: try: pref = reference_state['pressure_reference'] rgas = reference_state['ideal_gas_constant'] tref = reference_state['temperature_reference'] if self.gam is None: self.gamma = reference_state['specific_heat_ratio'].value except KeyError: vals = ('pressure_reference', 'ideal_gas_constant', 'temperature_reference', 'specific_heat_ratio') raise AttributeError('For temperature_stagnation, reference_state' ' is missing one or more of %s.' % (vals,)) if self.pressure is None: raise NotImplementedError('Get dimensional temperature_stagnation' ' from Q variables') rhoref = pref / rgas / tref vref = (rgas * tref).sqrt() momref = rhoref * vref self.pref = pref.value self.rgas = rgas.value self.rhoref = rhoref.value self.momref = momref.value # self.e0ref = ? self.tref = tref if cylindrical: self.mom_c1 = None if momentum.z is None else momentum.z.item self.mom_c2 = momentum.r.item self.mom_c3 = momentum.t.item else: self.mom_c1 = momentum.x.item self.mom_c2 = None if momentum.y is None else momentum.y.item self.mom_c3 = None if momentum.z is None else momentum.z.item
[docs] def calculate(self, loc, geom): """ Return metric value. """ rho = self.density(*loc) * self.rhoref vu = 0. if self.mom_c1 is None else self.mom_c1(*loc) * self.momref / rho vv = 0. if self.mom_c2 is None else self.mom_c2(*loc) * self.momref / rho vw = 0. if self.mom_c3 is None else self.mom_c3(*loc) * self.momref / rho if self.gam is not None: gamma = self.gam(*loc) else: gamma = self.gamma u2 = vu*vu + vv*vv + vw*vw if self.pressure is not None: ps = self.pressure(*loc) * self.pref else: e0 = self.energy(*loc) * self.e0ref / rho ps = (gamma-1.) * rho * (e0 - 0.5*u2) a2 = (gamma * ps) / rho mach2 = u2 / a2 ts = ps / (rho * self.rgas) return ts * (1. + (gamma-1.)/2. * mach2)
[docs] def dimensionalize(self, value): """ Dimensionalize `value`. """ return PhysicalQuantity(value, self.tref.get_unit_name())
register_metric('temperature_stagnation', TotalTemperature, False)
[docs]class Volume(object): """ Computes volume of mesh volume. """ def __init__(self, zone, zone_name, reference_state): if reference_state is None: self.units = None else: try: lref = reference_state['length_reference'] except KeyError: raise AttributeError("For volume, reference_state is missing" " 'length_reference'.") volref = lref * lref * lref self.units = volref.get_unit_name() self.volref = volref.value
[docs] def calculate(self, loc, volume): """ Return metric value. """ return volume * self.volref
[docs] def dimensionalize(self, value): """ Return dimensional `value`. """ return PhysicalQuantity(value, self.units)
register_metric('volume', Volume, True, 'volume')
OpenMDAO Home