jacobian_square_first Module Function

pure elemental module function jacobian_square_first(self, i, j, r) result(jacobian)

dx dx_dxi dx_deta dy dy_dxi dy_deta

Arguments

Type IntentOptional Attributes Name
class(type_square_first), intent(in) :: self
integer(kind=int32), intent(in) :: i
integer(kind=int32), intent(in) :: j
type(type_dp_vector_3d), intent(in) :: r

Return Value real(kind=real64)


Calls

proc~~jacobian_square_first~~CallsGraph proc~jacobian_square_first jacobian_square_first interface~dpsi_square_first type_square_first%dpsi_square_first proc~jacobian_square_first->interface~dpsi_square_first proc~get_coordinate abst_mesh%get_coordinate proc~jacobian_square_first->proc~get_coordinate proc~get_num_nodes~3 abst_mesh%get_num_nodes proc~jacobian_square_first->proc~get_num_nodes~3 proc~dpsi_square_first dpsi_square_first interface~dpsi_square_first->proc~dpsi_square_first

Called by

proc~~jacobian_square_first~~CalledByGraph proc~jacobian_square_first jacobian_square_first interface~jacobian_square_first type_square_first%jacobian_square_first interface~jacobian_square_first->proc~jacobian_square_first proc~is_in_square_first is_in_square_first proc~is_in_square_first->interface~jacobian_square_first interface~jacobian_det_square_first type_square_first%jacobian_det_square_first proc~is_in_square_first->interface~jacobian_det_square_first proc~jacobian_det_square_first jacobian_det_square_first proc~jacobian_det_square_first->interface~jacobian_square_first interface~is_in_square_first type_square_first%is_in_square_first interface~is_in_square_first->proc~is_in_square_first interface~jacobian_det_square_first->proc~jacobian_det_square_first proc~get_area_square_first get_area_square_first proc~get_area_square_first->interface~jacobian_det_square_first interface~get_area_square_first type_square_first%get_area_square_first interface~get_area_square_first->proc~get_area_square_first

Source Code

    pure elemental module function jacobian_square_first(self, i, j, r) result(jacobian)
        implicit none
        class(type_square_first), intent(in) :: self
        integer(int32), intent(in) :: i
        integer(int32), intent(in) :: j
        type(type_dp_vector_3d), intent(in) :: r
        real(real64) :: jacobian

        integer(int32) :: ii
        type(type_dp_vector_3d) :: coordinate

        jacobian = 0
        !! dx
        select case (i)
        case (1)
            select case (j)
            case (1)
                !! dx_dxi
                do ii = 1, self%get_num_nodes()
                    coordinate = self%get_coordinate(ii)
                    jacobian = jacobian + self%dpsi(ii, 1, r) * coordinate%x
                end do
            case (2)
                !! dx_deta
                do ii = 1, self%get_num_nodes()
                    coordinate = self%get_coordinate(ii)
                    jacobian = jacobian + self%dpsi(ii, 2, r) * coordinate%x
                end do
            end select

        !! dy
        case (2)
            select case (j)
            case (1)
                !! dy_dxi
                do ii = 1, self%get_num_nodes()
                    coordinate = self%get_coordinate(ii)
                    jacobian = jacobian + self%dpsi(ii, 1, r) * coordinate%y
                end do
            case (2)
                !! dy_deta
                do ii = 1, self%get_num_nodes()
                    coordinate = self%get_coordinate(ii)
                    jacobian = jacobian + self%dpsi(ii, 2, r) * coordinate%y
                end do
            end select
        end select

    end function jacobian_square_first