Limiter.py

# -*- coding: utf-8 -*-
"""
@author: Laura C. Kühle

"""
from abc import ABC, abstractmethod
import numpy as np


class Limiter(ABC):
    """Abstract class for limiting function.

    Methods
    -------
    get_name()
        Returns string of class name.
    apply(projection, cells)
        Applies limiting to cells.

    """
    def __init__(self, config):
        """Initializes Quadrature.

        Parameters
        ----------
        config : dict
            Additional parameters for limiter.

        """
        self._reset(config)

    @abstractmethod
    def _reset(self, config):
        """Resets instance variables.

        Parameters
        ----------
        config : dict
            Additional parameters for quadrature.

        """
        pass

    def get_name(self):
        """Returns string of class name."""
        return self.__class__.__name__

    @abstractmethod
    def apply(self, projection, cells):
        """Applies limiting to cells.

        Parameters
        ----------
        projection : ndarray
            Matrix of projection for each polynomial degree.
        cells : list
            Index of cells to limit.

        Returns
        -------
        ndarray
            Matrix of updated projection for each polynomial degree.

        """
        pass


class NoLimiter(Limiter):
    """Class without any limiting.

    Methods
    -------
    get_name()
        Returns string of class name.
    apply(projection, cells)
        Applies no limiting to cells.

    """
    def apply(self, projection, cells):
        """Returns projection without limiting."""
        return projection


class MinMod(Limiter):
    """Class for minmod limiting function.

    Sets projection for higher degrees to zero when forward backward, and cell
    slope values have not the same sign.

    Attributes
    ----------
    erase_degree : int
        Polynomial degree up to which projection is not set to zero during
        limiting.

    Methods
    -------
    get_name()
        Returns string of class name.
    apply(projection, cells)
        Applies limiting to cells.

    """
    def _reset(self, config):
        """Resets instance variables.

        Parameters
        ----------
        config : dict
            Additional parameters for quadrature.

        """
        # Unpack necessary configurations
        self._erase_degree = config.pop('erase_degree', 0)

    def get_name(self):
        """Returns string of class name concatenated with the erase-degree."""
        return self.__class__.__name__ + str(self._erase_degree)

    def apply(self, projection, cells):
        """Applies limiting to cells.

        Parameters
        ----------
        projection : ndarray
            Matrix of projection for each polynomial degree.
        cells : list
            Index of cells to limit.

        Returns
        -------
        new_projection : ndarray
            Matrix of updated projection for each polynomial degree.

        """
        new_projection = projection.copy()

        # If no troubled cells are detected, return copy
        if len(cells) == 0:
            return new_projection

        # Set mask to limit complete projection
        cell_slopes = self._set_cell_slope(projection)
        modification_mask = self._determine_mask(projection, cell_slopes)
        cells = np.array(cells)
        mask = np.zeros_like(new_projection, dtype=bool)
        mask[self._erase_degree+1:, cells+1] = np.tile(
            np.logical_not(modification_mask)[cells],
            (len(projection)-self._erase_degree-1, 1))

        # Limit troubled cells for higher degrees
        new_projection[mask] = 0

        return new_projection

    def _determine_mask(self, projection, slopes):
        """Determine limiting mask.

        Parameters
        ----------
        projection : ndarray
            Matrix of projection for each polynomial degree.
        slopes : ndarray
            Vector of slopes of projection cells.

        Returns
        -------
        ndarray
            Mask whether cells should be adjusted.

        """
        forward_slopes = (projection[0, 2:]-projection[0, 1:-1]) * (0.5**0.5)
        backward_slopes = (projection[0, 1:-1]-projection[0, :-2]) * (0.5**0.5)
        pos_mask = np.logical_and(slopes >= 0,
                                  np.logical_and(forward_slopes >= 0,
                                                 backward_slopes >= 0))
        neg_mask = np.logical_and(slopes <= 0,
                                  np.logical_and(forward_slopes <= 0,
                                                 backward_slopes <= 0))
        slope_mask = np.logical_or(pos_mask, neg_mask)

        return slope_mask

    @staticmethod
    def _set_cell_slope(projection):
        """Calculates the slope of the cell.

        Parameters
        ----------
        projection : ndarray
            Matrix of projection for each polynomial degree.

        Returns
        -------
        ndarray
            Vector of slopes of projection cells.

        """
        root_vector = np.array([np.sqrt(degree+0.5)
                                for degree in range(len(projection))])
        slope = root_vector[1:] @ projection[1:]
        return slope[1:-1]


class ModifiedMinMod(MinMod):
    """Class for modified minmod limiting function.

    Sets projection for higher degrees to zero when forward backward, and cell
    slope values have not the same sign and cell slope is significantly high.

    Attributes
    ----------
    threshold : float
        Threshold up to which a cell slope does not require limiting.

    Methods
    -------
    get_name()
        Returns string of class name.

    Notes
    -----
    Also called Cockburn-Shu limiter.

    """

    def _reset(self, config):
        """Resets instance variables.

        Parameters
        ----------
        config : dict
            Additional parameters for quadrature.

        """
        super()._reset(config)

        # Unpack necessary configurations
        cell_len = config.pop('cell_len')
        mod_factor = config.pop('mod_factor', 0)
        self._threshold = mod_factor * cell_len**2

    def get_name(self):
        """Returns string of class name concatenated with the erase-degree."""
        return self.__class__.__name__ + str(self._erase_degree)

    def _determine_mask(self, projection, slopes):
        """Determine limiting mask.

        Parameters
        ----------
        projection : ndarray
            Matrix of projection for each polynomial degree.
        slopes : ndarray
            Vector of slopes of projection cells.

        Returns
        -------
        ndarray
            Mask whether cells should be adjusted.

        """
        return np.logical_or(abs(slopes) <= self._threshold,
                             super()._determine_mask(projection, slopes))