2 #include "sll_working_precision.h" 
    3 #include "sll_assert.h" 
   14      sll_int32 :: n_dofs(2)
 
   17      sll_real64, 
allocatable :: conp(:,:)
 
   19      sll_real64 :: tl(3, 3)
 
   35     sll_int32, 
intent( in )  :: n_dofs(2) 
 
   36     sll_real64, 
intent( in ) :: conp(:,:) 
 
   37     sll_int32, 
intent( in ) :: k 
 
   38     sll_real64, 
intent( in ) :: tau 
 
   40     allocate( self%conp(product(n_dofs), 2) )
 
   47     self%tl(:,1) = 1._f64/3._f64
 
   48     self%tl(1,2) = 2._f64/(3._f64*self%tau)
 
   49     self%tl(2,2) = -1._f64/(3._f64*self%tau)
 
   50     self%tl(3,2) = -1._f64/(3._f64*self%tau)
 
   52     self%tl(2,3) = 1._f64/(sqrt(3._f64)*self%tau)
 
   53     self%tl(3,3) = -1._f64/(sqrt(3._f64)*self%tau)
 
   55     self%n_rows = product(self%n_dofs)
 
   56     self%n_cols = product(self%n_dofs) - self%n_dofs(2) + self%nk
 
   58     self%n_global_rows = self%n_rows
 
   59     self%n_global_cols = self%n_cols
 
   71     sll_real64, 
intent( in    ) :: x(:)
 
   72     sll_real64, 
intent(   out ) :: y(:)
 
   74     sll_int32 :: i, j, ind
 
   77     do j = 1, self%n_dofs(2)
 
   79           y(1+(j-1)*self%n_dofs(1)) = y(1+(j-1)*self%n_dofs(1)) + &
 
   80                ( self%conp(j, 1) * self%tl(i,2) + &
 
   81                self%conp(j, 2) * self%tl(i,3) )*&
 
   87     do j = 1, self%n_dofs(2)
 
   88        do i = 2, self%n_dofs(1)
 
   89           y(i+(j-1)*self%n_dofs(1)) =  x(ind)
 
   93     sll_assert(ind == self%n_cols+1)
 
module for abstract linear operator
class for abstract linear operator