Juan 07/03/2019: add posibilite of d/dx 2second order for Newmann , working...

Juan 07/03/2019: add posibilite of d/dx 2second order for Newmann , working well for 100 smoothing but not un global MG

tensorproductmultigrid_Source/discretisation.f90

@@ -293,32 +293,36 @@ contains
     integer, intent(in) :: iy
     real(kind=rl), intent(inout), dimension(5) :: alpha_T
     real(kind=rl), intent(out) :: Tij
+
+    !local var
     real(kind=rl)              :: h, rho_i, sigma_j
     real(kind=rl)              :: xcoef 
+    logical                    :: l2nd       
 
     h = grid_param%L/grid_param%n
     ! Cartesian coefficients
     Tij = h**2
     xcoef = 0.0
+    l2nd = .false.
+    alpha_T(1) = 1.0
+    alpha_T(2) = 1.0
     if (ix == grid_param%n) then
-      alpha_T(1) = 2.0_rl * xcoef
-    else
-      alpha_T(1) = 1.0_rl
+      alpha_T(1) = xcoef * 2.0_rl
+      if (l2nd) alpha_T(2) = 2.0_rl
     end if
     if (ix == 1) then
-      alpha_T(2) = 2.0_rl * xcoef
-    else
-      alpha_T(2) = 1.0_rl
+      alpha_T(2) = xcoef * 2.0_rl
+      if (l2nd) alpha_T(1) = 2.0_rl
     end if
+    alpha_T(3) = 1.0
+    alpha_T(4) = 1.0
     if (iy == grid_param%n) then
-      alpha_T(3) = 2.0_rl * xcoef
-    else
-      alpha_T(3) = 1.0_rl
+      alpha_T(3) = xcoef * 2.0_rl
+      if (l2nd) alpha_T(4) = 2.0
     end if
     if (iy == 1) then
-      alpha_T(4) = 2.0_rl * xcoef
-    else
-      alpha_T(4) = 1.0_rl
+      alpha_T(4) = xcoef * 2.0_rl
+      if (l2nd) alpha_T(3) = 2.0
     end if
 
     alpha_T(5) = alpha_T(1) + alpha_T(2) + alpha_T(3) + alpha_T(4)