"""Algorithms for building canonical data structure for integrals over subdomains."""
# Copyright (C) 2009-2016 Anders Logg and Martin Sandve Alnæs
#
# This file is part of UFL (https://www.fenicsproject.org)
#
# SPDX-License-Identifier: LGPL-3.0-or-later
import numbers
import typing
from collections import defaultdict
import ufl
from ufl.algorithms.coordinate_derivative_helpers import (
attach_coordinate_derivatives,
strip_coordinate_derivatives,
)
from ufl.form import Form
from ufl.integral import Integral
from ufl.protocols import id_or_none
from ufl.sorting import cmp_expr, sorted_expr
from ufl.utils.sorting import canonicalize_metadata, sorted_by_key
[docs]class IntegralData(object):
"""Utility class.
This class has members (domain, integral_type, subdomain_id,
integrals, metadata), where metadata is an empty dictionary that may
be used for associating metadata with each object.
"""
__slots__ = (
"domain",
"integral_type",
"subdomain_id",
"integrals",
"metadata",
"integral_coefficients",
"enabled_coefficients",
)
def __init__(self, domain, integral_type, subdomain_id, integrals, metadata):
"""Initialise."""
if 1 != len(set(itg.ufl_domain() for itg in integrals)):
raise ValueError("Multiple domains mismatch in integral data.")
if not all(integral_type == itg.integral_type() for itg in integrals):
raise ValueError("Integral type mismatch in integral data.")
if not all(subdomain_id == itg.subdomain_id() for itg in integrals):
raise ValueError("Subdomain id mismatch in integral data.")
self.domain = domain
self.integral_type = integral_type
self.subdomain_id = subdomain_id
self.integrals = integrals
# This is populated in preprocess using data not available at
# this stage:
self.integral_coefficients = None
self.enabled_coefficients = None
# TODO: I think we can get rid of this with some refactoring
# in ffc:
self.metadata = metadata
def __lt__(self, other):
"""Check if self is less than other."""
# To preserve behaviour of extract_integral_data:
return (self.integral_type, self.subdomain_id, self.integrals, self.metadata) < (
other.integral_type,
other.subdomain_id,
other.integrals,
other.metadata,
)
def __eq__(self, other):
"""Check for equality."""
# Currently only used for tests:
return (
self.integral_type == other.integral_type
and self.subdomain_id == other.subdomain_id
and self.integrals == other.integrals
and self.metadata == other.metadata
)
def __str__(self):
"""Format as a string."""
s = f"IntegralData over domain({self.integral_type}, {self.subdomain_id})"
s += " with integrals:\n"
s += "\n\n".join(map(str, self.integrals))
s += "\nand metadata:\n{metadata}"
return s
[docs]class ExprTupleKey(object):
"""Tuple comparison helper."""
__slots__ = ("x",)
def __init__(self, x):
"""Initialise."""
self.x = x
def __lt__(self, other):
"""Check if self is less than other."""
# Comparing expression first
c = cmp_expr(self.x[0], other.x[0])
if c < 0:
return True
elif c > 0:
return False
else:
# Comparing form compiler data
mds = canonicalize_metadata(self.x[1])
mdo = canonicalize_metadata(other.x[1])
return mds < mdo
[docs]def group_integrals_by_domain_and_type(integrals, domains):
"""Group integrals by domain and type.
Args:
integrals: list of Integral objects
domains: list of AbstractDomain objects from the parent Form
Returns:
Dictionary mapping (domain, integral_type) to list(Integral)
"""
integrals_by_domain_and_type = defaultdict(list)
for itg in integrals:
if itg.ufl_domain() is None:
raise ValueError("Integral has no domain.")
key = (itg.ufl_domain(), itg.integral_type())
# Append integral to list of integrals with shared key
integrals_by_domain_and_type[key].append(itg)
return integrals_by_domain_and_type
[docs]def integral_subdomain_ids(integral):
"""Get a tuple of integer subdomains or a valid string subdomain from integral."""
did = integral.subdomain_id()
if isinstance(did, numbers.Integral):
return (did,)
elif isinstance(did, tuple):
if not all(isinstance(d, numbers.Integral) for d in did):
raise ValueError("Expecting only integer subdomains in tuple.")
return did
elif did in ("everywhere", "otherwise"):
# TODO: Define list of valid strings somewhere more central
return did
else:
raise ValueError(f"Invalid domain id {did}.")
[docs]def rearrange_integrals_by_single_subdomains(
integrals: typing.List[Integral], do_append_everywhere_integrals: bool
) -> typing.Dict[int, typing.List[Integral]]:
"""Rearrange integrals over multiple subdomains to single subdomain integrals.
Args:
integrals: List of integrals
do_append_everywhere_integrals: Boolean indicating if integrals
defined on the whole domain should
just be restricted to the set of input subdomain ids.
Returns:
The integrals reconstructed with single subdomain_id
"""
# Split integrals into lists of everywhere and subdomain integrals
everywhere_integrals = []
subdomain_integrals = []
for itg in integrals:
dids = integral_subdomain_ids(itg)
if dids == "otherwise":
raise ValueError("'otherwise' integrals should never occur before preprocessing.")
elif dids == "everywhere":
everywhere_integrals.append(itg)
else:
subdomain_integrals.append((dids, itg))
# Fill single_subdomain_integrals with lists of integrals from
# subdomain_integrals, but split and restricted to single
# subdomain ids
single_subdomain_integrals = defaultdict(list)
for dids, itg in subdomain_integrals:
# Region or single subdomain id
for did in dids:
# Restrict integral to this subdomain!
single_subdomain_integrals[did].append(itg.reconstruct(subdomain_id=did))
# Add everywhere integrals to each single subdomain id integral
# list
otherwise_integrals = []
for ev_itg in everywhere_integrals:
# Restrict everywhere integral to 'otherwise'
otherwise_integrals.append(ev_itg.reconstruct(subdomain_id="otherwise"))
# Restrict everywhere integral to each subdomain
# and append to each integral list
if do_append_everywhere_integrals:
for subdomain_id in sorted(single_subdomain_integrals.keys()):
single_subdomain_integrals[subdomain_id].append(
ev_itg.reconstruct(subdomain_id=subdomain_id)
)
if otherwise_integrals:
single_subdomain_integrals["otherwise"] = otherwise_integrals
return single_subdomain_integrals
[docs]def build_integral_data(integrals):
"""Build integral data given a list of integrals.
The integrals you pass in here must have been rearranged and
gathered (removing the "everywhere" subdomain_id). To do this, you
should call group_form_integrals.
Args:
integrals: An iterable of Integral objects.
Returns:
A tuple of IntegralData objects.
"""
itgs = defaultdict(list)
# --- Merge integral data that has the same integrals,
unique_integrals = defaultdict(tuple)
metadata_table = defaultdict(dict)
for integral in integrals:
integrand = integral.integrand()
integral_type = integral.integral_type()
ufl_domain = integral.ufl_domain()
metadata = integral.metadata()
meta_hash = hash(canonicalize_metadata(metadata))
subdomain_id = integral.subdomain_id()
subdomain_data = id_or_none(integral.subdomain_data())
if subdomain_id == "everywhere":
raise ValueError(
"'everywhere' not a valid subdomain id. "
"Did you forget to call group_form_integrals?"
)
unique_integrals[(integral_type, ufl_domain, meta_hash, integrand, subdomain_data)] += (
subdomain_id,
)
metadata_table[(integral_type, ufl_domain, meta_hash, integrand, subdomain_data)] = metadata
for integral_data, subdomain_ids in unique_integrals.items():
(integral_type, ufl_domain, metadata, integrand, subdomain_data) = integral_data
integral = Integral(
integrand,
integral_type,
ufl_domain,
subdomain_ids,
metadata_table[integral_data],
subdomain_data,
)
# Group for integral data (One integral data object for all
# integrals with same domain, itype, (but possibly different metadata).
itgs[(ufl_domain, integral_type, subdomain_ids)].append(integral)
# Build list with canonical ordering, iteration over dicts
# is not deterministic across python versions
def keyfunc(item):
(d, itype, sid), integrals = item
sid_int = tuple(-1 if i == "otherwise" else i for i in sid)
return (d._ufl_sort_key_(), itype, (type(sid).__name__,), sid_int)
integral_datas = []
for (d, itype, sid), integrals in sorted(itgs.items(), key=keyfunc):
integral_datas.append(IntegralData(d, itype, sid, integrals, {}))
return integral_datas