See the following link

http://www.ijaamm.com/uploads/2/1/4/8/21481830/v2n1p3.pdf

Dear Dr. Bewar Ahamed,

      First I would like to thank for your answer. The paper you had mentioned having the contents of finding the exact solutions by means of analytical method. My point of interest is how to numerically solve NNV equations. Could you have any reference for that?

Thank You

The mentioned equations represent 2-dimensional generalization of (KdV) Korteweg-de-Vries equation u_t+uu_x+u_xxx=0, where subscript denotes partial derivative; see for example https://en.wikipedia.org/wiki/Korteweg%E2%80%93de_Vries_equation . Such equation has a solution in a form of soliton, i.e. in the functional form of a wave u(x-Vt), but with a unique peak. Solutions of similar structure were obtained in the paper recommended by Bewar.

As for numerical solution, one needs to start from  Cauchi problem with some initial condition. The knowledge of analytical solution allows to have an idea how solutions look like. But one can get a very large variety of behavior depending on initial conditions. Another problem is that such soliton propagates in unbounded space and one needs to select a proper boundary condition to avoid "reflection".

And here are some ideas for numerical solution of KdV: http://www.wikiwaves.org/Numerical_Solution_of_the_KdV . However, long ago I was using a different scheme for soliton, direct differences without Fourier transform. One also needs to select a good scheme so that there is no accumulation of numerical error that can destroy a solution. Similar method might work for the listed equation as well.

Hello Sudharsan Balasubramanian:        Your question is quite interesting and its numerical solution will be found in the same way as for the simple linear heat PDE and the more complicated nonlinear  Navier Stokes PDEs. The CRANK -NICOLSON time marching scheme as given by Galerkin Finite element discretization at an elemental level is 

                               D*Δn(t+δt)    =  E    .....................                 (1)

Where   using Matlab matrix symbolism                                                            D is the 3 by 3 matrix  D =    [K1 + 2P/δt,   K2,          K3 ;      Q,    -R,    0;     R ,      0 ,     -Q ]                                                                                              

 and E is the 3 by 1  column vector defined by E = [(2P/δt – K1)*Un(t) -  K2*Vn (t) – K3*Wn(t),     R* Vn(t) – Q*Un(t),   Q* Wn(t) – R*Un(t)  ]T,                                      Δn  = [Un,  Vn ,  Wn  ],                                                                               where Un, Vn and Wn are the nodal values of u, v, w ;                                 Un (t + δt)   =  nodal values of u at time = t + δt,  Un(t) = nodal values of u at time = t  ;                                                                                                      the matrices K1, K2, K3 are defined by  K1 = a*F + 3a*vav*Q + 3b*wav*R,               K2 = b*G + 3a*uav*Q,                        K3 =   3b*uav*R ;                                          uav = α*Un  ;  vav =   α*Vn,  wav = α*Wn,  and  α = (1/vol) ∫N dvol ;                     and the matrices P, Q, R, F, G are  defined by the integrals                                P =  ∫NT*N dvol,             Q =  ∫NT*Nx dvol,           R =  ∫NT*Ny dvol,                                F =   ∫NxT*Nxx dvol,  G =   ∫NyT*Nyy dvol,  

where Nx  is the 1st partial derivative of N wrt x.  In addition N is the shape function  row vector defined by N = [ N1  N2  N3  ...........................Nm] , m = number of nodes in an element,    

 with u = N*Un,  v = N*Vn,   w =  N*Wn,   Also Un = [Un1  Un2     Un3................Unm]  etc.  Un1 = value of u at node 1.              

  The discretization procedure can be found on page 37-39, and 86-88 of Smith and Griffith (1988) Programming the Finite Element Method. The CRANK-NICOLSON time marching scheme as given above may give an oscilatory solution. In that case it should be replaced by a FULLY IMPLICIT time marching scheme.

Hello Sudharsan Balasubramanian: 

   As mentioned the Crank-Nicholson time marching scheme may sometimes yield an oscillatory solution (op. cit. page 88). In that case a fully implicit scheme should be adopted.  This involves changes in the structure of the D and E matrices which shall now be as follows for the fully implicit scheme :

D is the 3 by 3 matrix  D =    [K1 + P/δt,       K2,          K3 ;      Q,        -R,         0;        R ,      0 ,        -Q ]                                                      

and E is the much simpler 3 by 1  column vector defined by E = [P*Un(t)/δt,       0,          0 ]T

where all the matrices and vectors retain their previous meanings.

The solution of the time marching equation shall be by iteration at each time step.  

Since, the average values uav, vav ,wav  are unknown at time t + δt  they have to be estimated before Eq. (1) of the previous post can be solved. A good choice is to assume Un (t + δt) ≈  Un (t ),  Vn (t + δt) ≈  Vn (t ), Wn (t + δt) ≈  Wn (t ),   compute the average values of u,v,w , substitute them i Eq. (1) and obtain new values of   Un (t + δt),  Vn (t + δt),  Wn (t + δt), compute the new average values of u, v, w, re-substitute into Eq.1 and solve Eq.(1) again. The iteration should stop when the nodal values of u, v, w stabilize and do not change further. A minimum error criterion should be used to spot when this occurs.

Dear Dr. Amaechi J. Anyaegbunam,

         Thanks for your answers. I'll go through it and back to you shortly. Once again I thank for your efforts.