!&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& 
SUBROUTINE lagrang3( x1, x2, x3, xst,  y1, y2, y3, yst) 
!***********************************************************************
!                                                                      *
!        Lagrangian 2nd order interpolation                            *
!                                                                      *
!***********************************************************************

      x1_x2 = x1 - x2
      x1_x3 = x1 - x3
      x2_x3 = x2 - x3

      xst_x1 = xst - x1
      xst_x2 = xst - x2
      xst_x3 = xst - x3

      p1 = xst_x1/x2_x3
      p2 = xst_x2/x1_x3
      p3 = xst_x3/x1_x2

      yst = p2 * p3 * y1  -  p1 * p3 * y2   +  p1 * p2 * y3

!-------------------------------------------------------------------------
END SUBROUTINE lagrang3