construct_square_second Module Function

module function construct_square_second(id, global_coordinate, input) result(element)

Arguments

Type IntentOptional Attributes Name
integer(kind=int32), intent(in) :: id
type(type_dp_3d), intent(in), pointer :: global_coordinate
type(type_input), intent(in) :: input

Return Value class(abst_element), allocatable


Calls

proc~~construct_square_second~~CallsGraph proc~construct_square_second construct_square_second interface~allocate_array allocate_array proc~construct_square_second->interface~allocate_array proc~initialize_abst_mesh abst_mesh%initialize_abst_mesh proc~construct_square_second->proc~initialize_abst_mesh proc~type_vtk_cell_get_dimension type_vtk_cell%type_vtk_cell_get_dimension proc~construct_square_second->proc~type_vtk_cell_get_dimension proc~type_vtk_cell_get_order type_vtk_cell%type_vtk_cell_get_order proc~construct_square_second->proc~type_vtk_cell_get_order proc~allocate_rank1_int16 allocate_rank1_int16 interface~allocate_array->proc~allocate_rank1_int16 proc~allocate_rank1_int32 allocate_rank1_int32 interface~allocate_array->proc~allocate_rank1_int32 proc~allocate_rank1_int64 allocate_rank1_int64 interface~allocate_array->proc~allocate_rank1_int64 proc~allocate_rank1_int8 allocate_rank1_int8 interface~allocate_array->proc~allocate_rank1_int8 proc~allocate_rank1_logical1 allocate_rank1_logical1 interface~allocate_array->proc~allocate_rank1_logical1 proc~allocate_rank1_logical4 allocate_rank1_logical4 interface~allocate_array->proc~allocate_rank1_logical4 proc~allocate_rank1_logical8 allocate_rank1_logical8 interface~allocate_array->proc~allocate_rank1_logical8 proc~allocate_rank1_real128 allocate_rank1_real128 interface~allocate_array->proc~allocate_rank1_real128 proc~allocate_rank1_real32 allocate_rank1_real32 interface~allocate_array->proc~allocate_rank1_real32 proc~allocate_rank1_real64 allocate_rank1_real64 interface~allocate_array->proc~allocate_rank1_real64 proc~allocate_rank2_int16 allocate_rank2_int16 interface~allocate_array->proc~allocate_rank2_int16 proc~allocate_rank2_int32 allocate_rank2_int32 interface~allocate_array->proc~allocate_rank2_int32 proc~allocate_rank2_int64 allocate_rank2_int64 interface~allocate_array->proc~allocate_rank2_int64 proc~allocate_rank2_int8 allocate_rank2_int8 interface~allocate_array->proc~allocate_rank2_int8 proc~allocate_rank2_logical1 allocate_rank2_logical1 interface~allocate_array->proc~allocate_rank2_logical1 proc~allocate_rank2_logical4 allocate_rank2_logical4 interface~allocate_array->proc~allocate_rank2_logical4 proc~allocate_rank2_logical8 allocate_rank2_logical8 interface~allocate_array->proc~allocate_rank2_logical8 proc~allocate_rank2_real128 allocate_rank2_real128 interface~allocate_array->proc~allocate_rank2_real128 proc~allocate_rank2_real32 allocate_rank2_real32 interface~allocate_array->proc~allocate_rank2_real32 proc~allocate_rank2_real64 allocate_rank2_real64 interface~allocate_array->proc~allocate_rank2_real64 proc~initialize_abst_mesh->interface~allocate_array none~set~3 type_dp_vector_3d%set proc~initialize_abst_mesh->none~set~3 proc~set_dp_vector_3d type_dp_vector_3d%set_dp_vector_3d none~set~3->proc~set_dp_vector_3d proc~set_dp_vector_3d_array type_dp_vector_3d%set_dp_vector_3d_array none~set~3->proc~set_dp_vector_3d_array proc~error_message error_message proc~allocate_rank1_int16->proc~error_message proc~allocate_rank1_int32->proc~error_message proc~allocate_rank1_int64->proc~error_message proc~allocate_rank1_int8->proc~error_message proc~allocate_rank1_logical1->proc~error_message proc~allocate_rank1_logical4->proc~error_message proc~allocate_rank1_logical8->proc~error_message proc~allocate_rank1_real128->proc~error_message proc~allocate_rank1_real32->proc~error_message proc~allocate_rank1_real64->proc~error_message proc~allocate_rank2_int16->proc~error_message proc~allocate_rank2_int32->proc~error_message proc~allocate_rank2_int64->proc~error_message proc~allocate_rank2_int8->proc~error_message proc~allocate_rank2_logical1->proc~error_message proc~allocate_rank2_logical4->proc~error_message proc~allocate_rank2_logical8->proc~error_message proc~allocate_rank2_real128->proc~error_message proc~allocate_rank2_real32->proc~error_message proc~allocate_rank2_real64->proc~error_message log_error log_error proc~error_message->log_error

Called by

proc~~construct_square_second~~CalledByGraph proc~construct_square_second construct_square_second interface~construct_square_second construct_square_second interface~construct_square_second->proc~construct_square_second interface~type_square_second type_square_second interface~type_square_second->interface~construct_square_second

Source Code

    module function construct_square_second(id, global_coordinate, input) result(element)
        implicit none
        integer(int32), intent(in) :: id
        type(type_dp_3d), pointer, intent(in) :: global_coordinate
        type(type_input), intent(in) :: input
        class(abst_element), allocatable :: element

        integer(int32) :: i
        integer(int32) :: num_nodes, num_gauss
        real(real64), allocatable :: weight(:)
        real(real64), allocatable :: gauss(:, :)

        allocate (type_square_second :: element)

        num_nodes = input%geometry%vtk%cells(id)%num_nodes_in_cell

        select case (input%basic%geometry_settings%integration_type)
        case ("full")
            num_gauss = 9_int32
            call allocate_array(weight, num_gauss)
            call allocate_array(gauss, 3_int32, num_gauss)

            weight(:) = [25.0d0 / 81.0d0, 40.0d0 / 81.0d0, 25.0d0 / 81.0d0, 40.0d0 / 81.0d0, &
                         64.0d0 / 81.0d0, 40.0d0 / 81.0d0, 25.0d0 / 81.0d0, 40.0d0 / 81.0d0, &
                         25.0d0 / 81.0d0]

            gauss(:, 1) = [-sqrt(3.0d0 / 5.0d0), -sqrt(3.0d0 / 5.0d0), 0.0d0]
            gauss(:, 2) = [0.0d0, -sqrt(3.0d0 / 5.0d0), 0.0d0]
            gauss(:, 3) = [sqrt(3.0d0 / 5.0d0), -sqrt(3.0d0 / 5.0d0), 0.0d0]
            gauss(:, 4) = [-sqrt(3.0d0 / 5.0d0), 0.0d0, 0.0d0]
            gauss(:, 5) = [0.0d0, 0.0d0, 0.0d0]
            gauss(:, 6) = [sqrt(3.0d0 / 5.0d0), 0.0d0, 0.0d0]
            gauss(:, 7) = [-sqrt(3.0d0 / 5.0d0), sqrt(3.0d0 / 5.0d0), 0.0d0]
            gauss(:, 8) = [0.0d0, sqrt(3.0d0 / 5.0d0), 0.0d0]
            gauss(:, 9) = [sqrt(3.0d0 / 5.0d0), sqrt(3.0d0 / 5.0d0), 0.0d0]
        case ("reduced")
            num_gauss = 4_int32
            call allocate_array(weight, num_gauss)
            call allocate_array(gauss, 3_int32, num_gauss)

            weight(:) = [1.0d0, 1.0d0, 1.0d0, 1.0d0]
            gauss(:, 1) = [-sqrt(1.0d0 / 3.0d0), -sqrt(1.0d0 / 3.0d0), 0.0d0]
            gauss(:, 2) = [-sqrt(1.0d0 / 3.0d0), sqrt(1.0d0 / 3.0d0), 0.0d0]
            gauss(:, 3) = [sqrt(1.0d0 / 3.0d0), sqrt(1.0d0 / 3.0d0), 0.0d0]
            gauss(:, 4) = [sqrt(1.0d0 / 3.0d0), -sqrt(1.0d0 / 3.0d0), 0.0d0]
        case ("free")
            num_gauss = 4_int32
            call allocate_array(weight, num_gauss)
            call allocate_array(gauss, 3_int32, num_gauss)

            weight(:) = [1.0d0, 1.0d0, 1.0d0, 1.0d0]
            gauss(:, 1) = [-input%basic%geometry_settings%integration_points, -input%basic%geometry_settings%integration_points, 0.0d0]
            gauss(:, 2) = [-input%basic%geometry_settings%integration_points, input%basic%geometry_settings%integration_points, 0.0d0]
            gauss(:, 3) = [input%basic%geometry_settings%integration_points, input%basic%geometry_settings%integration_points, 0.0d0]
            gauss(:, 4) = [input%basic%geometry_settings%integration_points, -input%basic%geometry_settings%integration_points, 0.0d0]
        end select

        call element%initialize(id=id, &
                                type=input%geometry%vtk%cells(id)%cell_type, &
                                group=input%geometry%vtk%cells(id)%cell_entity_id, &
                                dimension=input%geometry%vtk%cells(id)%get_dimension(), &
                                order=input%geometry%vtk%cells(id)%get_order(), &
                                num_nodes=num_nodes, &
                                connectivity=input%geometry%vtk%cells(id)%connectivity(1:num_nodes), &
                                num_gauss=num_gauss, &
                                weight=weight, &
                                gauss=gauss, &
                                global_coordinate=global_coordinate)

    end function construct_square_second