[cfdplus]

Hello,
Before everything I need to apologize for not writing the equation in a suitable format, with a quick search I came to know that I have to use Latex but unfortunately I have no experience with Latex.
I am trying to solve the transient stokes equation using finite element method. I discretized the equation using triangular element and ended up with below equation.

Stokes equation:

∂u/∂t=ϑ∆u-∂P/ρ

velocity x-direction:
(LNN+dt*nu*DNDN)*(u_x at time step n+1)=LNN*(u_x at time step n)+dt/rho(transpose of NDNx)*(Pressure at time step n+1)

velocity y-direction:
(LNN+dt*nu*DNDN)*(u_y at time step n+1)=LNN*(u_y at time step n)+dt/rho(transpose of NDNy)*(Pressure at time step n+1)+g*LNN

pressure equation:

DNDN*(pressure at time step n+1)=-(rho/dt)*(transpose of NDNx* (u_x at time step n+1)+transpose NDNy*(u_y at time step n+1)) + (Transpose NDNx* gradient P in x direction)+(transpose NDNy*gradient P in y direction))

I used the dt=0.001, nu= 1e-3, density=1000 and gravity=10. The problem is I get very large value for pressure and I'm notsure whether the discretized equations are correct ot not. The initial value for ux=2, uy=0, pressure=0. Imagine the domain is between -1<x<1, -1<y<1. The boundary condition is ux=2 for x=-1, uy=0 everywhere and p=1 at x=1.



LNN: lumped matrix of mass matrix
dt: time
nu:kinematic viscosity
rho:density
g:gravity
N : finite element shape function
Nx: derivation of shape function with respect to x
Ny: derivation of shape function with respect to y
NN:Mass Matrix
NDNx=Matrix of NNx in finite element
NDNy=Matrix of NNy in finite element
DNDN: Laplacian Matrix =-(NxNx+NyNy)