from math import sin, cos, pi

class doublegyre2d(object):
    """The Double Gyre vector field."""

    def __init__(self, a=0.1, eps=0.25, omega=pi/5):
        self.A = a
        self.EPS = eps
        self.OMEGA = omega

    # Samples the Double Gyre vector field
    def sample(self, x, y, t):
        u = -self.A * pi * sin(pi*(self.EPS*sin(self.OMEGA*t)*x*x + (1 - 2 * self.EPS*sin(self.OMEGA*t))*x))*cos(pi*y)
        v = self.A * pi*(2 * self.EPS*sin(self.OMEGA*t)*x - 2 * self.EPS*sin(self.OMEGA*t) + 1)* cos(pi*(self.EPS*sin(self.OMEGA*t)*x*x + (1 - 2 * self.EPS*sin(self.OMEGA*t))*x))*sin(pi*y)
        return u, v;

    # Samples the x-partial of the Double Gyre vector field
    def sample_dx(self, x, y, t):
        u = -self.A * pi*pi * (2 * self.EPS*sin(self.OMEGA*t)*x - 2 * self.EPS*sin(self.OMEGA*t) + 1)*cos(pi*(self.EPS*sin(self.OMEGA*t)*x*x + (1 - 2 * self.EPS*sin(self.OMEGA*t))*x))*cos(pi*y)
        v = 2 * self.A*self.EPS*pi* sin(self.OMEGA*t)*cos(pi*(self.EPS*sin(self.OMEGA*t)*x*x + (1 - 2 * self.EPS*sin(self.OMEGA*t))*x))*sin(pi*y) - self.A * pi*pi * (2 * self.EPS*sin(self.OMEGA*t)*x - 2 * self.EPS*sin(self.OMEGA*t) + 1) * (2 * self.EPS*sin(self.OMEGA*t)*x - 2 * self.EPS*sin(self.OMEGA*t) + 1) * sin(pi*(self.EPS*sin(self.OMEGA*t)*x*x + (1 - 2 * self.EPS*sin(self.OMEGA*t))*x))*sin(pi*y)
        return u, v;

    # Samples the y-partial of the Double Gyre vector field
    def sample_dy(self, x, y, t):
        u = self.A * pi*pi * sin(pi*(self.EPS*sin(self.OMEGA*t)*x*x + (1 - 2 * self.EPS*sin(self.OMEGA*t))*x))*sin(pi*y)
        v = self.A * pi*pi * (2 * self.EPS*sin(self.OMEGA*t)*x - 2 * self.EPS*sin(self.OMEGA*t) + 1)*cos(pi*(self.EPS*sin(self.OMEGA*t)*x*x + (1 - 2 * self.EPS*sin(self.OMEGA*t))*x))*cos(pi*y)
        return u, v;

    # Samples the t-partial of the Double Gyre vector field
    def sample_dt(self, x, y, t):
        u = -self.A * pi*pi * (self.EPS*self.OMEGA*cos(self.OMEGA*t)*x*x - 2 * self.EPS*self.OMEGA*cos(self.OMEGA*t)*x)*cos(pi*(self.EPS*sin(self.OMEGA*t)*x*x + (1 - 2 * self.EPS*sin(self.OMEGA*t))*x))*cos(pi*y)
        v = self.A * pi *	(2 * self.EPS*self.OMEGA*cos(self.OMEGA*t)*x - 2 * self.EPS*self.OMEGA*cos(self.OMEGA*t))*cos(pi*(self.EPS*sin(self.OMEGA*t)*x*x + (1 - 2 * self.EPS*sin(self.OMEGA*t))*x))*sin(pi*y) - self.A * pi*pi * (2 * self.EPS*sin(self.OMEGA*t)*x - 2 * self.EPS*sin(self.OMEGA*t) + 1)*(self.EPS*self.OMEGA*cos(self.OMEGA*t)*x*x - 2 * self.EPS*self.OMEGA*cos(self.OMEGA*t)*x) * sin(pi*(self.EPS*sin(self.OMEGA*t)*x*x + (1 - 2 * self.EPS*sin(self.OMEGA*t))*x))*sin(pi*y)
        return u, v;