Adding a custom model

If Adding a custom boundary condition is insufficient for your needs, it is also relatively straightforward to implement a custom model. In this case, however, the Discretization class will need to be modified (or the appropriate methods needs to be replaced in an instance of the class). The methods that most likely need to be changed are linear_part(), nonlinear_part(), mass_matrix() and boundaries(), or functions called from any of these functions.

As an example, the 2D nondimensional incompressible Stokes equations of the form

\[\begin{split}\frac{1}{\mathrm{Re}} \nabla^2 u - \nabla p &= 0\\ \nabla \cdot u &= 0\end{split}\]

could be implemented as

def linear_part(self):
    Re = self.get_parameter('Reynolds Number', 1.0)
    return 1 / Re * (self.u_xx() + self.u_yy()
                     + self.v_xx() + self.v_yy()) \
        - (self.p_x() + self.p_y()) \
        + self.div()

def nonlinear_part(self, state):
    atom = numpy.zeros((self.nx, self.ny, self.nz, self.dof, self.dof, 3, 3, 3))
    return atom, atom

We could then add a moving lid to the top using

def boundaries(self, atom):
    boundary_conditions = BoundaryConditions(
        self.nx, self.ny, self.nz, self.dim, self.dof, self.x, self.y, self.z)

    lid_velocity = self.parameters['Lid Velocity']
    boundary_conditions.moving_lid_north(atom, lid_velocity)

    boundary_conditions.no_slip_east(atom)
    boundary_conditions.no_slip_west(atom)

    boundary_conditions.no_slip_south(atom)

    return boundary_conditions.get_forcing()

\[\]