sll_m_linear_operator_ecsim_ev.F90

Go to the documentation of this file.

module sll_m_linear_operator_ecsim_ev
#include "sll_working_precision.h"
  use sll_m_linear_operator_abstract
  
  use sll_m_spline_fem_utilities, only : &
       sll_s_spline_fem_mass_line
  
  implicit none
  
  public :: sll_t_linear_operator_ecsim_ev
  
  private
  
  type, extends(sll_t_linear_operator_abstract) :: sll_t_linear_operator_ecsim_ev
     
     sll_int32 :: degree 
     sll_real64 :: dt 
     sll_real64 :: dx 
     sll_real64, allocatable :: mass_line_0(:) 
     sll_real64, allocatable :: mass_line_1(:) 
     sll_real64, pointer :: m1(:,:), m2(:,:), m4(:,:) 
 
 
   contains
     procedure :: create => create_linear_operator_ecsim_ev
     procedure :: free => free_ecsim_ev
     procedure :: dot  => dot_mono_r2r_ecsim_ev
     procedure :: print_info => print_info_ecsim_ev
 
  end type sll_t_linear_operator_ecsim_ev
 
contains
 
  subroutine create_linear_operator_ecsim_ev( self, n_dof, degree, &
       mass_line_0,mass_line_1, m1, m2, m4)
    class(sll_t_linear_operator_ecsim_ev), intent( inout ) :: self
    sll_int32 :: n_dof 
    sll_int32 :: degree 
    sll_real64 :: mass_line_0(:) 
    sll_real64 :: mass_line_1(:) 
    sll_real64,target :: m1(:,:) 
    sll_real64,target :: m2(:,:) 
    sll_real64,target :: m4(:,:) 
 
 
    self%n_dof=n_dof
    self%degree=degree
    self%dt=0.0_f64
    self%m1=>m1
    self%m2=>m2
    self%m4=>m4
    allocate(self%mass_line_0(degree+1))
    allocate(self%mass_line_1(degree))
    self%mass_line_0=mass_line_0
    self%mass_line_1=mass_line_1
 
 
    self%n_rows = 2*n_dof
    self%n_cols = 2*n_dof
 
    self%n_global_rows = self%n_rows
    self%n_global_cols = self%n_cols
  end subroutine create_linear_operator_ecsim_ev
 
 
  subroutine free_ecsim_ev( self )
    class(sll_t_linear_operator_ecsim_ev), intent( inout ) :: self
    deallocate(self%mass_line_0)
    deallocate(self%mass_line_1)
    self%m1=>null()
    self%m2=>null()
    self%m4=>null()
 
  end subroutine free_ecsim_ev
 
  subroutine dot_mono_r2r_ecsim_ev( self, x, y )
    class(sll_t_linear_operator_ecsim_ev), intent( in ) :: self
    sll_real64, intent( in    ) :: x(:) 
    sll_real64, intent(   out ) :: y(:) 
    ! local variables
    sll_int32 :: row, column
 
    y=0.0_f64
    do row=1, self%n_dof
       ! Multiplication of the inputvector with entries of the mass and particle-mass matrices
       y(row) = (self%mass_line_1(1)+0.5_f64*self%dt*self%m1(1,row))*x(row)
       y(self%n_dof+row)= (self%mass_line_0(1)+0.5_f64*self%dt*self%m4(1,row))*x(self%n_dof+row)
 
       do column = 2, self%degree
          y(row) = y(row) +&
               self%mass_line_1(column) * &
               (x(modulo(row+column-2,self%n_dof)+1) +&
               x(modulo(row-column,self%n_dof)+1))+&
               0.5_f64*self%dt*self%m1(column,row) * &
               x(modulo(row+column-2,self%n_dof)+1)+& 
               0.5_f64*self%dt*self%m1(column,modulo(row-column,self%n_dof)+1) * &
               x(modulo(row-column,self%n_dof)+1) 
 
       end do
       do column = 2, self%degree+1
          y(self%n_dof+row) = y(self%n_dof+row) +&
               self%mass_line_0(column) * &
               (x(self%n_dof+modulo(row+column-2,self%n_dof)+1) +&
               x(self%n_dof+ modulo(row-column,self%n_dof)+1))+&
               0.5_f64*self%dt*self%m4(column,row)*&
               x(self%n_dof+modulo(row+column-2,self%n_dof)+1)+&
               0.5_f64*self%dt*self%m4(column,modulo(row-column,self%n_dof)+1)*&
               x(self%n_dof+modulo(row-column,self%n_dof)+1)
       end do
 
       do column = 1, 2*self%degree
          y(row) = y(row) - &             
               0.5_f64*self%dt*self%m2(7-column,modulo(column-4+row-1,self%n_dof)+1)*&
               x(self%n_dof+modulo(row+column-self%degree-2,self%n_dof)+1) 
 
          y(row+self%n_dof) = y(row+self%n_dof) + &
               0.5_f64*self%dt*self%m2(column,row) * &
               x(modulo(row+column-self%degree-1,self%n_dof)+1)
       end do
    end do
 
  end subroutine dot_mono_r2r_ecsim_ev
 
  subroutine print_info_ecsim_ev( self )
    class(sll_t_linear_operator_ecsim_ev), intent(in) :: self
  end subroutine print_info_ecsim_ev
 
end module sll_m_linear_operator_ecsim_ev