Report Typos and Errors    
Semi-Lagrangian Library
Modular library for kinetic and gyrokinetic simulations of plasmas in fusion energy devices.
sll_m_linear_operator_schur_eb_1d.F90
Go to the documentation of this file.
2 #include "sll_working_precision.h"
3 
5 
6  use sll_m_matrix_csr, only: &
8 
9  use sll_m_spline_fem_utilities, only: &
12 
13  implicit none
14 
16 
17  private
19 
20  type(sll_t_matrix_csr), pointer :: mass0
21  type(sll_t_matrix_csr), pointer :: mass1
22  sll_int32 :: n_dofs
23  sll_real64 :: delta_x
24  sll_real64 :: sign = 1.0_f64
25 
26  contains
27  procedure :: create => create_linear_operator_schur_eb_1d
28  procedure :: free => free_schur_eb_1d
29  procedure :: dot => dot_schur_eb_1d
30  procedure :: print_info => print_info_schur_eb_1d
31 
33 
34 
35 contains
36 
37  subroutine create_linear_operator_schur_eb_1d( self, mass0, mass1, n_dofs, delta_x )
38  class(sll_t_linear_operator_schur_eb_1d), intent( inout ) :: self
39  type(sll_t_matrix_csr), target :: mass0
40  type(sll_t_matrix_csr), target :: mass1
41  sll_int32 :: n_dofs
42  sll_real64 :: delta_x
43 
44  self%mass0 => mass0
45  self%mass1 => mass1
46  self%n_dofs = n_dofs
47  self%delta_x = delta_x
48 
49  self%n_rows = self%n_dofs
50  self%n_cols = self%n_dofs
51 
52  self%n_global_rows = self%n_rows
53  self%n_global_cols = self%n_cols
54 
56 
57  subroutine free_schur_eb_1d( self )
58  class(sll_t_linear_operator_schur_eb_1d), intent( inout ) :: self
59 
60  self%mass0=> null()
61  self%mass1=> null()
62 
63  end subroutine free_schur_eb_1d
64 
65 
66  subroutine dot_schur_eb_1d ( self, x, y )
67  class(sll_t_linear_operator_schur_eb_1d), intent( in ) :: self
68  sll_real64, intent( in ) :: x(:)
69  sll_real64, intent( out ) :: y(:)
70  !local variables
71  sll_real64 :: scratch1(self%n_dofs), scratch2(self%n_dofs)
72 
73  ! Compute D x
74  call sll_s_multiply_g_1d( self%n_dofs, self%delta_x, x, scratch1)
75 
76  ! Compute M1 D x
77  call self%mass1%dot( scratch1, scratch2 )
78 
79  ! Compute D^T M1 D x
80  call sll_s_multiply_gt_1d( self%n_dofs, self%delta_x, scratch2, scratch1 )
81 
82  ! Compute M0 x
83  call self%mass0%dot( x, y )
84 
85  ! Sum up the two parts
86  y = y + self%sign * scratch1
87 
88  end subroutine dot_schur_eb_1d
89 
90 
91  subroutine print_info_schur_eb_1d( self )
92  class(sll_t_linear_operator_schur_eb_1d), intent(in) :: self
93 
94  end subroutine print_info_schur_eb_1d
95 
96 
module for abstract linear operator
subroutine create_linear_operator_schur_eb_1d(self, mass0, mass1, n_dofs, delta_x)
module for Compressed Sparse Row Matrix (CSR)
Utilites for Maxwell solver's with spline finite elements.
subroutine, public sll_s_multiply_g_1d(n_dofs, delta_x, in, out)
Multiplication of the input vector in by the derivative matrix G.
subroutine, public sll_s_multiply_gt_1d(n_dofs, delta_x, in, out)
Multiplication of the input vector in by the transposed derivative matrix G^T
    Report Typos and Errors