mform/m_mform_constraint_abstract.f90 - GCC Code Coverage Report

Directory:	src/
Coverage:	low: ≥ 0% medium: ≥ 75.0% high: ≥ 90.0%

	Coverage	Exec / Excl / Total
Lines:	82.6%	38 / 0 / 46
Functions:	66.7%	2 / 0 / 3
Branches:	42.1%	53 / 0 / 126

    mform/m_mform_constraint_abstract.f90
    
        Line
        Branch
        Exec
        Source
      
        !> @brief Module containing an object representing a local linear constraint (i.e., one that acts on only part of the m-form) that
      
        !>        can be imposed on an m-form function (MFormFun)
      
        !>
      
        !> We consider a homogeneous linear constraint of the form
      
        !>   C * u[inds] = 0,
      
        !> is considered, where u is the vector of unknowns of a tensor product B-spline, C is a wide matrix and inds is a list of indices
      
        !> on which the constraint is applied.
      
        !>
      
        !> Two methods for imposing the constraint are considered:
      
        !> - using an extraction matrix E which spans the nullspace of C:
      
        !>    C * E = 0.
      
        !> - using a projector P which projects the coefficients onto the nullspace of C:
      
        !>   P = E (F^T * E)^(-1) F^T,
      
        !> where F is a matrix for which F^T * E is invertible.
      
        !>
      
        !> The current implementation restricts the domain decomposition in such a way that the local constraint acts on the coefficients of
      
        !> an m-form that are local to the current rank
      
        module m_mform_constraint_abstract
      
          use m_bspline, only: BSplineFun
      
          use m_common, only: wp
      
          use m_mform_basis, only: MFormSpace
      
          use m_sparsemat, only: SparseMat
      
          use m_tensorprod_indices, only: TensorProdIndexList
      
          use m_coord_transform_abstract, only: CoordTransformAbstract
      
        #include "petsc.fi"
      
          implicit none
      
          private
      
          public :: MFormConstraintLocal
      
          !> The threshold under which extraction/projection coefficients are considered zero
      
          real(wp), parameter :: COEFF_THR = 1e-14_wp
      
          !> @brief The object that represents linear constrains imposed locally (i.e., only acting on inds) on the coefficients of an mform
      
          type, abstract :: MFormConstraintLocal
      
            !> The m-form space to which the constraint applies
      
            type(MFormSpace) :: space
      
            !> Whether the constraint is homogeneous (currently always .true.)
      
            logical :: is_homogeneous
      
            !> Whether the projector is symmetric
      
            logical :: is_symmetric
      
            !> The indices of the coefficients that the constraint acts upon (per m-form component)
      
            type(TensorProdIndexList), allocatable :: index_list(:)
      
            !> Block-wise (per m-form component) sparse representation of the projector (restricted to the coefficients that the constraint
      
            !> acts upon)
      
            type(SparseMat), allocatable :: projector(:, :)
      
            !> Block-wise (per m-form component) sparse representation of the extraction operator (restricted to the coefficients that the
      
            !> constraint acts upon)
      
            type(SparseMat), allocatable :: extractor(:)
      
            !> If the constraint '.not. is_homogeneous', the prescribed values of the constraint are stored here per component
      
            ! type(VectorWrapper), allocatable :: inhomogeneous_values(:)
      
            !> Whether the constraint is shared across multiple MPI ranks (per spatial direction)
      
            logical :: is_mpi_shared(3)
      
          contains
      
            procedure(mform_constraint_local_init), deferred :: init
      
            procedure(mform_constraint_local_destroy), deferred :: destroy_derived
      
            procedure :: destroy => destroy_mform_constraint_local
      
            procedure :: apply => apply_mform_constraint
      
          end type
      
          interface
      
            subroutine mform_constraint_local_init(this, space, coord_transform)
      
              import :: MFormConstraintLocal, MFormSpace, CoordTransformAbstract
      
              class(MFormConstraintLocal), intent(inout) :: this
      
              type(MFormSpace), intent(in) :: space
      
              class(CoordTransformAbstract), intent(in), optional :: coord_transform
      
            end subroutine mform_constraint_local_init
      
            subroutine mform_constraint_local_destroy(this)
      
              import :: MFormConstraintLocal
      
              class(MFormConstraintLocal), intent(inout) :: this
      
            end subroutine mform_constraint_local_destroy
      
          end interface
      
        contains
      
          !> @brief Apply a constraint to an m-form
      
          !>
      
          !> @param[in] this The constraint to apply
      
          !> @param[out] vout The MForm to which the constraint has been applied
      
          !> @param[in] vin The MForm to which the constraint is applied
      
          4/4✓ Branch 2 → 3 taken 42 times.
✓ Branch 2 → 4 taken 2 times.
✓ Branch 4 → 5 taken 42 times.
✓ Branch 4 → 6 taken 2 times.

        44
          subroutine apply_mform_constraint(this, vout, vin)
      
            use m_mform_basis, only: MFormFun
      
            implicit none
      
            class(MFormConstraintLocal), intent(in) :: this
      
            type(MFormFun), intent(out) :: vout
      
            type(MFormFun), intent(in) :: vin
      
            integer :: blk_out, blk_in, local_idx, local_row, local_col
      
            integer :: i_in, j_in, k_in, i_out, j_out, k_out, row_start, row_end
      
            real(wp) :: local_val
      
          1/2✗ Branch 6 → 7 not taken.
✓ Branch 6 → 8 taken 44 times.

        44
            if (this%space%m /= vin%space%m) then
      
        ✗
              error stop "MFormConstraintLocal::apply_mform_constraint: MForm mismatch"
      
            end if
      
          1/2✗ Branch 8 → 9 not taken.
✓ Branch 8 → 10 taken 44 times.

        44
            if (.not. allocated(this%projector)) then
      
        ✗
              error stop "MFormConstraintLocal::apply_mform_constraint: projector not allocated"
      
            end if
      
            ! Make a deep copy
      
        44
            vout = vin
      
          2/2✓ Branch 12 → 13 taken 44 times.
✓ Branch 12 → 15 taken 92 times.

        136
            do blk_out = 1, size(vout%tp_funs)
      
          2/2✓ Branch 16 → 17 taken 24678 times.
✓ Branch 16 → 18 taken 92 times.

        24814
              do local_idx = 0, size(this%index_list(blk_out)%lin) - 1
      
        24678
                i_out = this%index_list(blk_out)%i(local_idx)
      
        24678
                j_out = this%index_list(blk_out)%j(local_idx)
      
        24678
                k_out = this%index_list(blk_out)%k(local_idx)
      
        24770
                vout%tp_funs(blk_out)%data(i_out, j_out, k_out) = 0._wp
      
              end do
      
            end do
      
          2/2✓ Branch 20 → 21 taken 92 times.
✓ Branch 20 → 34 taken 44 times.

        136
            do blk_out = 1, size(vout%tp_funs)
      
          2/2✓ Branch 22 → 23 taken 236 times.
✓ Branch 22 → 32 taken 92 times.

        372
            do blk_in = 1, size(vin%tp_funs)
      
          2/2✓ Branch 24 → 25 taken 61494 times.
✓ Branch 24 → 30 taken 236 times.

        61822
              do local_row = 0, this%projector(blk_out, blk_in)%nr_rows - 1
      
        61494
                row_start = this%projector(blk_out, blk_in)%row_starts_at_nz(local_row)
      
        61494
                row_end = this%projector(blk_out, blk_in)%row_starts_at_nz(local_row + 1) - 1
      
        61494
                i_out = this%index_list(blk_out)%i(local_row)
      
        61494
                j_out = this%index_list(blk_out)%j(local_row)
      
        61494
                k_out = this%index_list(blk_out)%k(local_row)
      
          1/2✓ Branch 26 → 27 taken 61494 times.
✗ Branch 26 → 29 not taken.

        61730
                do local_idx = row_start, row_end
      
        ✗
                  local_col = this%projector(blk_out, blk_in)%col_idx(local_idx)
      
        ✗
                  local_val = this%projector(blk_out, blk_in)%values(local_idx)
      
        ✗
                  i_in = this%index_list(blk_in)%i(local_col)
      
        ✗
                  j_in = this%index_list(blk_in)%j(local_col)
      
        ✗
                  k_in = this%index_list(blk_in)%k(local_col)
      
                  vout%tp_funs(blk_out)%data(i_out, j_out, k_out) = vout%tp_funs(blk_out)%data(i_out, j_out, k_out) &
      
        61494
                                                                    + local_val * vin%tp_funs(blk_in)%data(i_in, j_in, k_in)
      
                end do
      
              end do
      
            end do
      
            end do
      
        44
            call vout%distribute()
      
        44
          end subroutine apply_mform_constraint
      
          !> @brief Destroy the abstract MFormConstraintLocal object
      
          !>
      
          !> @param[inout] this The MFormConstraintLocal object to destroy
      
        2388
          subroutine destroy_mform_constraint_local(this)
      
            implicit none
      
            class(MFormConstraintLocal), intent(inout) :: this
      
            integer :: s1, s2
      
          2/2✓ Branch 2 → 3 taken 908 times.
✓ Branch 2 → 23 taken 1480 times.

        2388
            if (allocated(this%index_list)) then
      
          2/2✓ Branch 4 → 5 taken 1210 times.
✓ Branch 4 → 7 taken 908 times.

        2118
              do s1 = 1, size(this%index_list)
      
        2118
                call this%index_list(s1)%destroy()
      
              end do
      
          8/14✓ Branch 8 → 9 taken 908 times.
✗ Branch 8 → 20 not taken.
✓ Branch 10 → 11 taken 1210 times.
✓ Branch 10 → 20 taken 908 times.
✗ Branch 11 → 12 not taken.
✓ Branch 11 → 13 taken 1210 times.
✗ Branch 13 → 14 not taken.
✓ Branch 13 → 15 taken 1210 times.
✗ Branch 15 → 16 not taken.
✓ Branch 15 → 17 taken 1210 times.
✗ Branch 17 → 18 not taken.
✓ Branch 17 → 19 taken 1210 times.
✗ Branch 20 → 21 not taken.
✓ Branch 20 → 22 taken 908 times.

        2118
              deallocate (this%index_list)
      
            end if
      
          2/2✓ Branch 23 → 24 taken 908 times.
✓ Branch 23 → 42 taken 1480 times.

        2388
            if (allocated(this%extractor)) then
      
          2/2✓ Branch 25 → 26 taken 1210 times.
✓ Branch 25 → 28 taken 908 times.

        2118
              do s1 = 1, size(this%extractor)
      
        2118
                call this%extractor(s1)%destroy()
      
              end do
      
          7/12✓ Branch 29 → 30 taken 908 times.
✗ Branch 29 → 39 not taken.
✓ Branch 31 → 32 taken 1210 times.
✓ Branch 31 → 39 taken 908 times.
✗ Branch 32 → 33 not taken.
✓ Branch 32 → 34 taken 1210 times.
✗ Branch 34 → 35 not taken.
✓ Branch 34 → 36 taken 1210 times.
✗ Branch 36 → 37 not taken.
✓ Branch 36 → 38 taken 1210 times.
✗ Branch 39 → 40 not taken.
✓ Branch 39 → 41 taken 908 times.

        2118
              deallocate (this%extractor)
      
            end if
      
          2/2✓ Branch 42 → 43 taken 331 times.
✓ Branch 42 → 65 taken 2057 times.

        2388
            if (allocated(this%projector)) then
      
          2/2✓ Branch 44 → 45 taken 433 times.
✓ Branch 44 → 51 taken 331 times.

        764
              do s1 = 1, size(this%projector, 1)
      
          2/2✓ Branch 46 → 47 taken 739 times.
✓ Branch 46 → 49 taken 433 times.

        1503
                do s2 = 1, size(this%projector, 2)
      
        1172
                  call this%projector(s1, s2)%destroy()
      
                end do
      
              end do
      
          7/12✓ Branch 52 → 53 taken 331 times.
✗ Branch 52 → 62 not taken.
✓ Branch 54 → 55 taken 739 times.
✓ Branch 54 → 62 taken 331 times.
✗ Branch 55 → 56 not taken.
✓ Branch 55 → 57 taken 739 times.
✗ Branch 57 → 58 not taken.
✓ Branch 57 → 59 taken 739 times.
✗ Branch 59 → 60 not taken.
✓ Branch 59 → 61 taken 739 times.
✗ Branch 62 → 63 not taken.
✓ Branch 62 → 64 taken 331 times.

        1070
              deallocate (this%projector)
      
            end if
      
        2388
            call this%destroy_derived()
      
        2388
          end subroutine destroy_mform_constraint_local
      
        ✗
        end module m_mform_constraint_abstract

Line	Branch	Exec	Source
1			!> @brief Module containing an object representing a local linear constraint (i.e., one that acts on only part of the m-form) that
2			!> can be imposed on an m-form function (MFormFun)
3			!>
4			!> We consider a homogeneous linear constraint of the form
5			!> C * u[inds] = 0,
6			!> is considered, where u is the vector of unknowns of a tensor product B-spline, C is a wide matrix and inds is a list of indices
7			!> on which the constraint is applied.
8			!>
9			!> Two methods for imposing the constraint are considered:
10			!> - using an extraction matrix E which spans the nullspace of C:
11			!> C * E = 0.
12			!> - using a projector P which projects the coefficients onto the nullspace of C:
13			!> P = E (F^T * E)^(-1) F^T,
14			!> where F is a matrix for which F^T * E is invertible.
15			!>
16			!> The current implementation restricts the domain decomposition in such a way that the local constraint acts on the coefficients of
17			!> an m-form that are local to the current rank
18			module m_mform_constraint_abstract
19			use m_bspline, only: BSplineFun
20			use m_common, only: wp
21			use m_mform_basis, only: MFormSpace
22			use m_sparsemat, only: SparseMat
23			use m_tensorprod_indices, only: TensorProdIndexList
24			use m_coord_transform_abstract, only: CoordTransformAbstract
25			#include "petsc.fi"
26
27			implicit none
28
29			private
30			public :: MFormConstraintLocal
31
32			!> The threshold under which extraction/projection coefficients are considered zero
33			real(wp), parameter :: COEFF_THR = 1e-14_wp
34
35			!> @brief The object that represents linear constrains imposed locally (i.e., only acting on inds) on the coefficients of an mform
36			type, abstract :: MFormConstraintLocal
37			!> The m-form space to which the constraint applies
38			type(MFormSpace) :: space
39
40			!> Whether the constraint is homogeneous (currently always .true.)
41			logical :: is_homogeneous
42
43			!> Whether the projector is symmetric
44			logical :: is_symmetric
45
46			!> The indices of the coefficients that the constraint acts upon (per m-form component)
47			type(TensorProdIndexList), allocatable :: index_list(:)
48
49			!> Block-wise (per m-form component) sparse representation of the projector (restricted to the coefficients that the constraint
50			!> acts upon)
51			type(SparseMat), allocatable :: projector(:, :)
52
53			!> Block-wise (per m-form component) sparse representation of the extraction operator (restricted to the coefficients that the
54			!> constraint acts upon)
55			type(SparseMat), allocatable :: extractor(:)
56
57			!> If the constraint '.not. is_homogeneous', the prescribed values of the constraint are stored here per component
58			! type(VectorWrapper), allocatable :: inhomogeneous_values(:)
59
60			!> Whether the constraint is shared across multiple MPI ranks (per spatial direction)
61			logical :: is_mpi_shared(3)
62			contains
63			procedure(mform_constraint_local_init), deferred :: init
64			procedure(mform_constraint_local_destroy), deferred :: destroy_derived
65			procedure :: destroy => destroy_mform_constraint_local
66			procedure :: apply => apply_mform_constraint
67			end type
68
69			interface
70			subroutine mform_constraint_local_init(this, space, coord_transform)
71			import :: MFormConstraintLocal, MFormSpace, CoordTransformAbstract
72			class(MFormConstraintLocal), intent(inout) :: this
73			type(MFormSpace), intent(in) :: space
74			class(CoordTransformAbstract), intent(in), optional :: coord_transform
75			end subroutine mform_constraint_local_init
76
77			subroutine mform_constraint_local_destroy(this)
78			import :: MFormConstraintLocal
79			class(MFormConstraintLocal), intent(inout) :: this
80			end subroutine mform_constraint_local_destroy
81			end interface
82
83			contains
84
85			!> @brief Apply a constraint to an m-form
86			!>
87			!> @param[in] this The constraint to apply
88			!> @param[out] vout The MForm to which the constraint has been applied
89			!> @param[in] vin The MForm to which the constraint is applied
90	4/4 ✓ Branch 2 → 3 taken 42 times. ✓ Branch 2 → 4 taken 2 times. ✓ Branch 4 → 5 taken 42 times. ✓ Branch 4 → 6 taken 2 times.	44	subroutine apply_mform_constraint(this, vout, vin)
91			use m_mform_basis, only: MFormFun
92
93			implicit none
94
95			class(MFormConstraintLocal), intent(in) :: this
96			type(MFormFun), intent(out) :: vout
97			type(MFormFun), intent(in) :: vin
98
99			integer :: blk_out, blk_in, local_idx, local_row, local_col
100			integer :: i_in, j_in, k_in, i_out, j_out, k_out, row_start, row_end
101			real(wp) :: local_val
102
103	1/2 ✗ Branch 6 → 7 not taken. ✓ Branch 6 → 8 taken 44 times.	44	if (this%space%m /= vin%space%m) then
104		✗	error stop "MFormConstraintLocal::apply_mform_constraint: MForm mismatch"
105			end if
106
107	1/2 ✗ Branch 8 → 9 not taken. ✓ Branch 8 → 10 taken 44 times.	44	if (.not. allocated(this%projector)) then
108		✗	error stop "MFormConstraintLocal::apply_mform_constraint: projector not allocated"
109			end if
110
111			! Make a deep copy
112		44	vout = vin
113
114	2/2 ✓ Branch 12 → 13 taken 44 times. ✓ Branch 12 → 15 taken 92 times.	136	do blk_out = 1, size(vout%tp_funs)
115	2/2 ✓ Branch 16 → 17 taken 24678 times. ✓ Branch 16 → 18 taken 92 times.	24814	do local_idx = 0, size(this%index_list(blk_out)%lin) - 1
116		24678	i_out = this%index_list(blk_out)%i(local_idx)
117		24678	j_out = this%index_list(blk_out)%j(local_idx)
118		24678	k_out = this%index_list(blk_out)%k(local_idx)
119
120		24770	vout%tp_funs(blk_out)%data(i_out, j_out, k_out) = 0._wp
121			end do
122			end do
123
124	2/2 ✓ Branch 20 → 21 taken 92 times. ✓ Branch 20 → 34 taken 44 times.	136	do blk_out = 1, size(vout%tp_funs)
125	2/2 ✓ Branch 22 → 23 taken 236 times. ✓ Branch 22 → 32 taken 92 times.	372	do blk_in = 1, size(vin%tp_funs)
126	2/2 ✓ Branch 24 → 25 taken 61494 times. ✓ Branch 24 → 30 taken 236 times.	61822	do local_row = 0, this%projector(blk_out, blk_in)%nr_rows - 1
127		61494	row_start = this%projector(blk_out, blk_in)%row_starts_at_nz(local_row)
128		61494	row_end = this%projector(blk_out, blk_in)%row_starts_at_nz(local_row + 1) - 1
129
130		61494	i_out = this%index_list(blk_out)%i(local_row)
131		61494	j_out = this%index_list(blk_out)%j(local_row)
132		61494	k_out = this%index_list(blk_out)%k(local_row)
133
134	1/2 ✓ Branch 26 → 27 taken 61494 times. ✗ Branch 26 → 29 not taken.	61730	do local_idx = row_start, row_end
135		✗	local_col = this%projector(blk_out, blk_in)%col_idx(local_idx)
136		✗	local_val = this%projector(blk_out, blk_in)%values(local_idx)
137
138		✗	i_in = this%index_list(blk_in)%i(local_col)
139		✗	j_in = this%index_list(blk_in)%j(local_col)
140		✗	k_in = this%index_list(blk_in)%k(local_col)
141
142			vout%tp_funs(blk_out)%data(i_out, j_out, k_out) = vout%tp_funs(blk_out)%data(i_out, j_out, k_out) &
143		61494	+ local_val * vin%tp_funs(blk_in)%data(i_in, j_in, k_in)
144			end do
145			end do
146			end do
147			end do
148
149		44	call vout%distribute()
150		44	end subroutine apply_mform_constraint
151
152			!> @brief Destroy the abstract MFormConstraintLocal object
153			!>
154			!> @param[inout] this The MFormConstraintLocal object to destroy
155		2388	subroutine destroy_mform_constraint_local(this)
156			implicit none
157
158			class(MFormConstraintLocal), intent(inout) :: this
159
160			integer :: s1, s2
161
162	2/2 ✓ Branch 2 → 3 taken 908 times. ✓ Branch 2 → 23 taken 1480 times.	2388	if (allocated(this%index_list)) then
163	2/2 ✓ Branch 4 → 5 taken 1210 times. ✓ Branch 4 → 7 taken 908 times.	2118	do s1 = 1, size(this%index_list)
164		2118	call this%index_list(s1)%destroy()
165			end do
166	8/14 ✓ Branch 8 → 9 taken 908 times. ✗ Branch 8 → 20 not taken. ✓ Branch 10 → 11 taken 1210 times. ✓ Branch 10 → 20 taken 908 times. ✗ Branch 11 → 12 not taken. ✓ Branch 11 → 13 taken 1210 times. ✗ Branch 13 → 14 not taken. ✓ Branch 13 → 15 taken 1210 times. ✗ Branch 15 → 16 not taken. ✓ Branch 15 → 17 taken 1210 times. ✗ Branch 17 → 18 not taken. ✓ Branch 17 → 19 taken 1210 times. ✗ Branch 20 → 21 not taken. ✓ Branch 20 → 22 taken 908 times.	2118	deallocate (this%index_list)
167			end if
168	2/2 ✓ Branch 23 → 24 taken 908 times. ✓ Branch 23 → 42 taken 1480 times.	2388	if (allocated(this%extractor)) then
169	2/2 ✓ Branch 25 → 26 taken 1210 times. ✓ Branch 25 → 28 taken 908 times.	2118	do s1 = 1, size(this%extractor)
170		2118	call this%extractor(s1)%destroy()
171			end do
172	7/12 ✓ Branch 29 → 30 taken 908 times. ✗ Branch 29 → 39 not taken. ✓ Branch 31 → 32 taken 1210 times. ✓ Branch 31 → 39 taken 908 times. ✗ Branch 32 → 33 not taken. ✓ Branch 32 → 34 taken 1210 times. ✗ Branch 34 → 35 not taken. ✓ Branch 34 → 36 taken 1210 times. ✗ Branch 36 → 37 not taken. ✓ Branch 36 → 38 taken 1210 times. ✗ Branch 39 → 40 not taken. ✓ Branch 39 → 41 taken 908 times.	2118	deallocate (this%extractor)
173			end if
174	2/2 ✓ Branch 42 → 43 taken 331 times. ✓ Branch 42 → 65 taken 2057 times.	2388	if (allocated(this%projector)) then
175	2/2 ✓ Branch 44 → 45 taken 433 times. ✓ Branch 44 → 51 taken 331 times.	764	do s1 = 1, size(this%projector, 1)
176	2/2 ✓ Branch 46 → 47 taken 739 times. ✓ Branch 46 → 49 taken 433 times.	1503	do s2 = 1, size(this%projector, 2)
177		1172	call this%projector(s1, s2)%destroy()
178			end do
179			end do
180	7/12 ✓ Branch 52 → 53 taken 331 times. ✗ Branch 52 → 62 not taken. ✓ Branch 54 → 55 taken 739 times. ✓ Branch 54 → 62 taken 331 times. ✗ Branch 55 → 56 not taken. ✓ Branch 55 → 57 taken 739 times. ✗ Branch 57 → 58 not taken. ✓ Branch 57 → 59 taken 739 times. ✗ Branch 59 → 60 not taken. ✓ Branch 59 → 61 taken 739 times. ✗ Branch 62 → 63 not taken. ✓ Branch 62 → 64 taken 331 times.	1070	deallocate (this%projector)
181			end if
182
183		2388	call this%destroy_derived()
184
185		2388	end subroutine destroy_mform_constraint_local
186		✗	end module m_mform_constraint_abstract
187