Skip to content

Vectorized Bruss #34

@YingboMa

Description

@YingboMa
 kernel_u! = let N=N, xyd=xyd_brusselator, dx=step(xyd_brusselator)
    @inline function (du, u, A, B, α, II, I, t)
      i, j = Tuple(I)
      x = xyd[I[1]]
      y = xyd[I[2]]
      ip1 = limit(i+1, N); im1 = limit(i-1, N)
      jp1 = limit(j+1, N); jm1 = limit(j-1, N)
      du[II[i,j,1]] = α[II[i,j,1]]*(u[II[im1,j,1]] + u[II[ip1,j,1]] + u[II[i,jp1,1]] + u[II[i,jm1,1]] - 4u[II[i,j,1]])/dx^2 +
      B[II[i,j,1]] + u[II[i,j,1]]^2*u[II[i,j,2]] - (A[II[i,j,1]] + 1)*u[II[i,j,1]] + brusselator_f(x, y, t)
    end
  end
  kernel_v! = let N=N, xyd=xyd_brusselator, dx=step(xyd_brusselator)
    @inline function (du, u, A, B, α, II, I, t)
      i, j = Tuple(I)
      ip1 = limit(i+1, N)
      im1 = limit(i-1, N)
      jp1 = limit(j+1, N)
      jm1 = limit(j-1, N)
      du[II[i,j,2]] = α[II[i,j,2]]*(u[II[im1,j,2]] + u[II[ip1,j,2]] + u[II[i,jp1,2]] + u[II[i,jm1,2]] - 4u[II[i,j,2]])/dx^2 +
      A[II[i,j,1]]*u[II[i,j,1]] - u[II[i,j,1]]^2*u[II[i,j,2]]
    end
  end
  brusselator_2d = let N=N, xyd=xyd_brusselator, dx=step(xyd_brusselator)
    function (du, u, p, t)
      @inbounds begin
        ii1 = N^2
        ii2 = ii1+N^2
        ii3 = ii2+2(N^2)
        A = @view p[1:ii1]
        B = @view p[ii1+1:ii2]
        α = @view p[ii2+1:ii3]
        II = LinearIndices((N, N, 2))
        kernel_u!.(Ref(du), Ref(u), Ref(A), Ref(B), Ref(α), Ref(II), CartesianIndices((N, N)), t)
        kernel_v!.(Ref(du), Ref(u), Ref(A), Ref(B), Ref(α), Ref(II), CartesianIndices((N, N)), t)
        return nothing
      end
    end
  end

Metadata

Metadata

Assignees

Labels

No labels
No labels

Type

No type

Projects

No projects

Milestone

No milestone

Relationships

None yet

Development

No branches or pull requests

Issue actions